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Introduction 


Estimation of the sensitivity of problem functions with respect to problem variables 
forms the basis for many of our modem day algorithms for engineering optimization. The 
most common application of problem sensitivities has been in the calculation of objective 
function and constraint partial derivatives for determining se arch directions and optimality 
conditions. A second form of sensitivity analysis, paramete r sensitivity, has also become 
an important topic in recent years with the advent of renewed research in the optimization of 
large engineering systems by means of decomposition methods. By parameter sensitivity, 
we refer to the estimation of changes in the modeling functions and current design variables 
due to small changes in the fixed parameters of the formula! ion. Methods for calculating 
these derivatives have been proposed and have been used a ; the basis of a method for 
multi-level decomposition of large engineering problems [Mobieski, 1982], Two 
drawbacks to estimating parameter sensitivities by current methods have been. (1) the need 
for second order information about the Lagrangian at the current point, and (2) the 
estimates assume no change in the active set of constraints. The objectives of this work 
were to investigate solutions to these two problems. 


1 1 STANDARD NOTATION 


To provide a framework about which we can discuss the various ways sensitivity 
analysis can be performed, the following standard form of the nonlinear programming 
problem, which explicidy represents the problem parameters, is presented. 


Minimize: f(x,P) 

Subject to: hi(x,P) = 0 
gj(x,P) > 0 

x min - x - X max 
x = (xi,X2,...,x n ) 
P = (pi,P2.»'Pk) 


Objective function 
Equality constraints 1 = Id- 
Inequality constrtints j = 1,J 
Variable bounds 
Design variables 
Problem parameu rs 


( 1 . 1 ) 

( 1 . 2 ) 

(1.3) 

(1.4) 

(1.5) 

(1.6) 


In the above formulation, we assume that the problem functions f, g, and h can be 
either linear or nonlinear functions of the design variables. We also assume that the 
problem parameters P, are held fixed during the course of the optimization. Any candidate 
solution point, x* , must satisfy the following first order <uhn-Tucker conditions: 
V x L(x,v,u) = 0 
hi(x)=0 1 = 1 ,L 

gj(x)>0 j = U t (L9) 



ujgj(x)=0 j = l,J 

uj > 0 j = 1»J 

where the Lagrangian L, is given by: 

L(x,v,u) = f(x) + Svi hj(x) - luj gj(x) 

At some point, usually the optimal point, we are interested in understanding the 
effect that changes in P will have on our proposed solution >*. Therefore we seek the 
sensitivities, df/dP, dx/3P, and aOhgVdP 1 . In this report, ve will propose a new 
algorithm based on the Recursive Quadratic Programming (F QP) method for estimating 
these parameter sensitivities. The following sections provide a description of this algorithm 
and how it relates to current methods, a discussion of the implementation issues, and some 
initial testing on a test set of known characteristics. In addit on, section 6 proposes some 
solutions for estimating sensitivities in those cases where the active set of the constraints 
changes when the parameter is changed. 


( 1 . 10 ) 

( 1 . 11 ) 

( 1 . 12 ) 


1 The notation (h,g) refers to the set of constraints active at the curren point 
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2 . Background 


The standard problem of parameter sensitivity analysts is to indicate how the 
objective function, constraints, and optimum design variables will change when problem 
parameters or design variables are changed from their current values. Parameter Sensitivity 
analysis is usually performed at a candidate optimum point where we might be interested in 
studying how the optimal design might be effected by changes in specifications, variability 
due to manufacturing, or operational noises. In this chapter we present a historical 
overview of the significant developments in sensitivity analysis and provide a review and 
assessment of current parameter sensitivity methods. The final section of the chapter 
reviews work done in estimating parameter sensitivities for hose cases where the active 

constraint set changes. 

2,1, RF.VTEW OF PARAMETER SENSrTTV TTY METT ODS 

The roots of sensitivity analysis can be traced to Lagrange (1881) when he ^ 
suggested solving equality constrained extrema problems by finding the solution x*,and 

v*, for the equations 

V x L(x,v) = 0 (2 ' 1} 


where 

L(x,v) = f(x) + ]Lvihi(x) 

where the v t are undetermined multipliers or Lagrange multipliers. The paper did not 
provide the conditions for when solutions of equation (2. 1-2.3), were actual solutions of 
the extrema problems or how to interpret the Lagrange multipliers. 

Samuelson (1947) gave several interpretations of l^agrange multipliers in an 
economic setting. He developed approaches based on using Lagrange multipliers to solve 
different economic models and was the first to clearly identify Lagrange multipliers as 
shadow prices in an economic context. Kuhn and Tucker (1951) presented conditions for 
relative extrema which use the Lagrange multipliers to establish optimality (ref. eq. 1.7 - 
1.12). Since 1951 several constraint qualifications and extensions to these conditions have 
been proposed and are described in Bazaraa and Shetty (1 979). 

Dantzig (1963) brought forth the idea of "Post Optimality Analysis" for linear 
programs. Dantzig described post optimality analysis as the calculation of the sensitivity of 
the optimum with respect to changes in the problem parameters. Sensitivity analysis has 
been widely used in linear programming, a good survey of its use is provided by Gal 



(1984). 


Fiacco et al. (1968,1974,1976,1983) has also done e: . tensive research in the area of 
sensitivity analysis. His book "Introduction to Sensitivity ar d Stability Analysis" (1983) 
covers the significant developments in the field of sensitivit) analysis pnor to 1982. He 
has published many articles on sensitivity analysis, and has probably been the most active 
researcher of sensitivity analysis for nonlinear programming problems. 

In the following subsections, we will discuss past w< »rk related to the determination 
of sensitivity information for nonlinear programming problems. The methods we will 
discuss range from the most simplistic approach of reoptimisation to more elaborate 
approaches based on the Kuhn-Tucker conditions or advanc ed optimization methods. 

9, . 1 . 1 . Brute Force Methods 

The simplest, and probably most used method, for parameter sensitivity analysis is 
to re-optimize the problem for the new values of the problem parameters and plot the 
trends. We will refer to this as the Brute Force method. Ti e Brute Force method is 
probably the most accurate of the methods available (for large variations in Ap, but can 
experience round off and truncation errors when used to approximate derivatives) but it can 
be computationally expensive even for small problems. Ex imples of its use in the literature 
are given in Arbuckle and Sliwa (1984) and Robertson and Gabriele (1987). 

Armacost and Fiacco (1974) and McKeown (1980 b) describe a direct approach to 
calculating parameter sensitivities based on the central difference approximation given 
below 


df* 

f(x* ,p + Ap) - f(x*,p - Ap) 

(2.4) 

dp 

2Ap 

dx* 

x*(p + Ap) - x*(p - Ap) 

(2.5) 

w = 

2Ap 


This method requires the problem to be reoptimized (to a high degree of accuracy) for two 
different values of the parameter. McKeown states that th s method should not be used as a 
primary method for the calculation of sensitivities because it is computationally expensive. 

9 1 2 Knhn-Tucker Methods 

To avoid the computational expense of reoptimiza ion, several researchers have 
developed sensitivity methods based on the Kuhn-Tuckei conditions (1.7) - (1.12). Two 
types of algorithms have resulted, those that differentiate he Kuhn-Tucker conditions with 
respect to p, and those that differentiate the optimality cor ditions for penalty functions. 
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In the former category, a set of Kuhn-Tucker sensitiv ity equations have been 
derived independently by several authors (Armacost and Fia sco 1974, Sobieski et. at. 
1981, McKeown 1980 b) and result in the following linear system of equations. 


— o 

dx 


~dV x L~ 

L -V x (h,g) 

dpi 

+ 

dpi 

. V x (h,g)T 0 . 

d(v,u) 

- dpi - 


d(h,g) 

- dpi _ 


This linear system can be solved for the sensitivity c f the design variables with 
respect to a problem parameter dx/dpi, and the sensitivity of the Lagrange multipliers with 
respect to Pi , d(v,u)/dpi. These can then be used to determi tie the sensitivity of the 
objective function with respect to pi by the following 

df__ df df T dx (2.7) 

dpi dpi 3x dpi" 

For any change in the parameter A Pi , the new optimum value of the objective 
function or design variables can be estimated from the linear extrapolations 

,, * , A . df (2.8) 

fnew — f( x old) + ^Pi ^p. 

* * A dx (2.9) 

X new = * old + APi 


These equations are bounded by the assumption that the active set remains the 
same. An estimate of when the active set will change can 1 « made by examining the 
Lagrange multipliers of the active inequality constraints and linear approximations of the 
inactive constraints. An inequality constraint should leave the active set when its Lagrange 
multiplier goes to zero. The corresponding value of A Pi v here this occurs is predicted by 

using the linear prediction 


Api = 



lapiJ 


j e active set of constraints 


( 2 . 10 ) 


A new inequality constraint will enter the active set when its value goes to zero. A linear 
prediction for when this happens is given by 


a JL j 4 active set of cons raints 

(dgj dgj T dx^j 

Wi + dx 


( 2 . 11 ) 
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We can predict the change in active set to occur at the smalle .t value of Apt obtained from 
applying equations 2.10 and 2.1 1 to all constraints. 

Fiacco (1974,1980,1983) has developed first and set ond order extrapolation 
techniques to predict the new value of the optimum when parameters are perturbed. 
Armacost and Fiacco have developed a second order extrapolation for the objective function 
value for the special case where the problem parameters are confined to being the right hand 
side values of the constraints. This provides second order response information for the 
objective function using the Lagrange multipliers and the pa tials with respect to P of the 

Lagrange multipliers. 

Sobieski, et. al. (1981) observed that a more accurate estimate of W g lven m 
(2.8) can be obtained if the value of x n ew given in (2.9) is u ed to calculate the value of the 
objective function at a perturbation A Pi . This will be a mon accurate estimate for problems 
where the constraints are well behaved and not highly nonli tear, but the objective function 
is nonlinear. 

Barthelemy and Sobieski (1983) derived the follow ng formula that can also be 
used to calculate the sensitivity of the objective function without the need to calculate 

3x*/3p, 

nineq 

< 212 > 

dp op i op 

j =1 


The formula can be derived by assuming that objective function behaves like the 
Lagrangian in the region of the optimum. This formula ha s also been derived by Fiacco 

(1983) and McKeown (1980 b). 

Diewart (1984) has developed some new sensitivity theories for dealing with the 
addition of constraints at the solution of economic models before the solution of the 
sensitivity equations. This analysis is important because there may be short term 
restrictions on modifications that can be made to the system. The paper presents a 
recursive relationship that can be used to avoid refactorizi ig the sensitivity equations when 
a new constraint is added to the problem. The paper also presents equations that can be 
used to calculate a second order estimate of the location o ' the optimum, but this formula 
requires third order derivatives which are seldom available in engineering. 

7 13 Methods Based on the Exten ded Desien Sr ac£ 

Vanderplaats (1984 a, 1984 b) and Vanderplaats tnd Yoshida (1985, 1986) have 
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developed an approach for calculating the sensitivity based oi the method of feasible 
directions. The sensitivities are estimated by extending the sc t of design variables to 
include the problem parameters for which a feasible direction is then determined. This 
method is known as the Extended Design Space (EDS) method. Of the methods discussed, 
it has the dual advantages of simplicity and efficiency. Vanderplaats (1984 a) reports that 
the EDS method can handle near active constraints, and is ab e to leave constraint 
linearizations. However, the method does suffer from a sensitivity to one of its algorithm 
parameters as reported in Vanderplaats and Cai (1987), and is unable to predict when 
constraints will leave the active set. The EDS method is also sensitive to the restriction of 

the move vector to be of length one. 

The EDS method can be used to assess the effect of 3 lerturbing several parameters at 
the same time. It is also able to solve for sensitivities of degenerate optimal points where 
either strict complementarity does not hold, or the constraint gradients of the active 
constraints are linearly dependent. The method seems to gh e good estimations for medium 
sized perturbations of the parameters, but for small perturbations the the Kuhn-Tucker 
method described above gives better results. Vanderplaats and Cai (1987) also report that 
there are some cases where the EDS algorithm can produce incorrect values of the 

sensitivity derivatives. 

Vanderplaats also proposes a second order approximation technique which is 
interesting but requires second derivatives of the objective function and constraints. The 
second order method solves a quadratic approximating problem for a specified value of the 
parameter. The second order method will give good result . in a larger region about the 
optimum than first order methods and does not appear to as sensitive to changes in the 
active set as other methods are. However, there is still the problem of obtaining the 
Hessians of the objective function and constraints and solving the quadratic approximating 
problem. Vanderplaats and Cai (1987) feel that the second order EDS algorithm is the best 
option short of reoptimizing the problem for estimating sensitivities. But they caution that 
the method should not always be used because of its high computational cost. 

2 } 4 Variable Sensitivities 

McKeown (1980 a,c) has developed sensitivity an llysis techniques for determining 
the sensitivity of design variables subject to perturbations about the optimum. This 
technique is based on an eigenvector analysis of the reduced Hessian matrix which applies 
to a variant of our standard problem (1.1)-(1.6) where nc problem parameters exist . For 
unconstrained problems the major eigenvector will point in the direction of maximum 
increase of the objective function, the minor eigenvector will point in the direction of 
minimum increase of the objective function. For constrained problems the directions are 



projected on the active constraints. This type of information may be useful for setting 
tolerances on design variables. 

For McKeown's algorithm, the Hessian of the Lagrangian is needed but the 
analysis is performed using only the reduced Hessian of the Lagrangian. An algorithm is 
provided for reducing the Hessian. If the Hessian is to be evaluated numerically, an 
algorithm is provided for the calculation of the reduced Hessian of the Lagrangian directly. 
This wiU reduce the number of extra function evaluations that axe needed to conduct the 

sensitivity analysis. 

7 15. Other Work 

Garcia and Zangwill (1981) describe a Homotopy approach that can be used to 
solve nonlinear programming problems. They state that this approach can also be used to 
solve parametric nonlinear programming problems and is closely related to sensitivity and 
perturbation analysis. Komija and Hirabay (1984) discuss some theoretical topics involved 
in using a Homotopy approach to calculate parameter sensi tivities when the active set of 

constraints changes. 

Dinkel and Kochenberger and Wong (1983) have developed an incremental 
approach for solving for the sensitivities of geometric prog ramming problems. The 
approach is to ask the user for the new value of the parameter and then make several steps 
with corrections to reach that point. They found the smaller the step they used the more 
accurate the solution would be. 

Jittomtrum (1984) examines solving for the sensiti vity of degenerate optimum 
points using the Kuhn-Tucker sensitivity equations. He provides a way to solve these 
problems using directional derivatives which provides different answers for both positive 
and negative perturbations in the parameters. Other theoretical issues for the use of 
directional derivatives to calculate optimum parameter sensitivities have been addressed by 
Janin (1984), Gauvin and Dubeau (1983), and Rockafell ir, R. T. (1984). 

Zolezzi (1985) examines the conditions under which the Lagrange multipliers are 
continuous under perturbations in the problem data. This is important because Kuhn- 
Tucker sensitivity analysis uses Lagrange multipliers and rates of change of the Lagrange 
multipliers to predict the rate of change of the objective function. Comet and Laroque 
(1987) establish conditions under which the values of the Lagrange multipliers are 
Lipschitz continuous for perturbations in the problem dat i. 

Ganesh and Biegler (1987) have developed a sensitivity analysis based on the 
reduced Hessian. The reduction is conducted by using tl e equality constraints and the 
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implicit function theorem to reduce the dimensionality of the Hessian matrix that needs to 
be calculated. Their method is beneficial when there are eqi ality constraints present in the 
formulation of the problem, because they have reduced the i lumber of function evaluations 
required to find the required second order information numc rically. Their method does not 
provide dv/9p without calculating the fuU Hessian of the Lapangian. 

Rao (1987 a) and Guang-Yaun and Wen-Quan (1985) have studied the problem of 
dealing with fuzzy constraints and fuzzy objective function ;. In their work they first solve 
a crisp problem then they attempt to calculate how far they < an relax constraints while 
improving the objective function. To use their technique the user is required to specify 
how much violation is allowed in the constraints. Templenian (1987) reports using fuzzy 
set theory and optimization to design structures and deal wi :h uncertainties in the problem. 

Sandgren, Gim and Ragsdell (1985) describe a pro! lem formulation that can be 
used to obtain optimum designs with a minimum sensitivity to uncontrollable parameters. 
Their approach does not use post optimality analysis but uses a modified objective function 
to deal with the uncertainties in the problem parameters* 

The area of calculating sensitivity derivatives with espect to design variables ( i. e. 
the calculation of gradients of functions) has been an area of active research. This can lead 
to significant savings over using finite differencing. The structural optimization community 
now widely uses sensitivity analysis when the finite element method is used to analyze a 
structure. An excellent survey article of methods of sensitivity analysis for structural 
optimization is provided by Adelman and Haftka (1986). 

Haug and Arora, et al. (1977,1979,1981) have de' eloped ways to calculate the 
gradients analytically for many structural and dynamic applications. Many of these 
methods are described in the book by Haug, Komkov and Choi (1985). 

Sobieski, et al. (1981,1982,1983,1984,1985,1986,1987) has been working on 
developing sensitivity techniques for use with multi-level decomposition techniques. 
Decomposition methods break the solution of a large prob em into a system level problem 
and a group of subproblems. Each subproblem is solved using a special formulation and 
inputs from the system level problem. A sensitivity analysis is performed on the 
subproblem and the results are feed as input to the system level problem. The system level 
problem gathers all the sensitivities of the subproblems and then based on these inputs and 
others, determines the next iteration of the process. Usually, the equations (2.6) - (2.9) are 
used at the subsystem level to determine the required sensitivities, but some difficulties 
have been encountered when changes in the active set occur. 
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Schmit and Chang (1984) have developed an extension of Sobieski's work and 
derived sensitivity equations for structural optimization problems. They derived more 
restrictive limits on the allowable perturbations than those provided by Sobieski. They 
have assumed that second derivatives of the constraints are av ailable which is true of many 
structural problems but may not be true for other application ;ireas. 

Schmit and Chang formulated their structural optimization problem using reciprocal 
variables and solved for the sensitivity of the dual problem. For their structural problems, 
the Hessian of the Lagrangian was diagonally dominate and t he Hessian of the objective 
function was analytically available. For this class of problems good results can be expected 
even if the Hessian of the Lagrangian is inaccurate. 

Buys and Gonin (1977) developed and implemented a sensitivity analysis 
procedure for an augmented Lagrangian (AL) type code, VF01A. Their implementation is 
encouraging because they make use of the approximations o‘ the Hessian of the Lagrangian 
that were calculated during the solution of the original problem, The results that they 
obtained using the approximate matrices were in very close agreement of those obtained by 

using the exact matrices. 

McKeown (1980 b) derives both the first and second order Kuhn-Tucker parameter 
sensitivity equations. He also provides a discussion of Fia< co's sensitivity for SUMT 
penalty functions versus Buys and Gonin’s sensitivity for A L penalty functions. He 
concludes that using sensitivity for AL penalty functions should be superior to sensitivity 
by SUMT because AL produces better conditioned matrice: . 

7 2 PRFVTOT IS WORK IN ESTIM A TING PARAMF FR SF.NSmvmES FQ R 
CHANGES IN THE ACTIVE $ E1 

When the active set of constraints changes, one of t he underlying assumptions 
made in deriving the Kuhn-Tucker sensitivity equations is violated. This can result in 
inaccuracies in any extrapolations based on these sensitivities since, in general, a change in 
the active constraints will result in a different set of sensitivities. Accurate sensitivity 
analysis in the presence of active set changes is also very ii nportant for efficient 
convergence of the multi-level decomposition techniques ] .reposed by Sobieski and, in 
general, for an accurate representation of the local sensitiv ties. 

In the following subsections, we will first discuss he different cases that occur as a 
result of a constraint entering or leaving the active set. wh.,t effects these cases have on 
sensitivity analysis, and how changes in the active set can be predicted. We will then 
present examples of the sensitivities for the different cases which will also serve to indicate 

how the different sensitivity algorithms perform. 
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2 2 - - 1 - Cases to Consider 


When a new constraint enters the active set, or a curr ently active constraint leaves 
the active set, we can expect a change in the sensitivity derivatives. However, it is also 
possible that the linear independence of the constraint gradit nts can also be affected. For 
the discussion that follows, we define the following four cates that can result from changes 


in the active set. 


1. A constraint enters the active set and the constraint gradients are linearly 
^Aconstraint leaves the active set and the constraint gradients are linearly 


^Aconsttaint enters the active set replacing an act ve constraint and the constraint 

gradients are linearly dependent. .. 

4. A constraint enters the active set and feasible region disappears. 


For Cases 1 and 2, we can expect discontinuities in the following derivatives when 
the active set changes: d2f*/dp2, dx*/dp, and 5u*/5p. 

Case 3 is characterized by a discontinuity in the Laf range multiplier estimates which 
causes a discontinuity in df*/dp. Since the active set changes there will also be a 
discontinuity in dx*/dp. At the point where the constraints become linearly dependent , the 
Kuhn-Tucker sensitivity equations become singular. Often what is happening for Case 3 is 
that an exchange of constraints in the active set is about tc take place (i.e. the new 
constraint may replace one of the constraints that is already in the active set ). If the 
problem is not poorly formulated, we will find ourselves moving through the degenerate 
point as p increases or decreases and one of the constraint will be dropped from the active 

set. 

Case 4 is characterized as a point from which p ca i only be perturbed in one 
direction. If p is perturbed in the wrong direction this wil cause there to be no feasible 
region and there will be no solution for the optimization problem with this value of p. Thus 
we can only perturb p in the one direction that causes the optimum path to move into the 
feasible region, and there will only exist a directional derivative for the problem in that 
direction. Case 4 can be thought of as an overconstraineo design where the designer 
adjusted a parameter to the point where the design is no longer able to meet specifications. 

7 2 2. Prediction of when the Activ e set will Cht.ng£ 

Barthelemy and Sobieski (1983 a) have observed that the accuracy of extrapolations 
of the objective function deteriorates rapidly when the ac ive set changes. From section 
2.1.2, we saw that we can use equations 2.10 and 2.1 1 to predict where the active set will 
change, thus we can use this information to predict when the extrapolations will deteriorate. 
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A problem with bounding Ap by equations 2. 10 and 2.1 1 is that the estimate is only 
good for the first constraint that is encountered because once the active set changes the 
search direction to the new optimum will change (the discoi tinuity in dx/dp). Thus, it 
becomes very difficult to estimate when or which constraint will leave/enter the active set 
second. This problem will be addressed in section 6. 


The merit of using equations 2. 10 and 2.1 1 to predit t when the active set will 
change was discussed by Adelman and Haftka (1986). They state, "The effectiveness of 
using this approach (equations 2.10 and 2.1 1) is still in douot with positive results being 
obtained by Schmit and Chang (1984) and negative results lieing obtained by Barthelemy 
and Sobieski (1983 a)". We feel that the positive results thi t were obtained by Schmit and 
Chang are due to problem linearity and the changes in the at tive set that they encountered 
being case 1 and case 2 changes. We feel that the negative -esults obtained by Barthelemy 
and Sobieski are due to nonlinearity of the problem and also a case 3 change in the active 
set taking place. As we will see later in this report, the con >equences on sensitivity 
derivatives of case 3 changes in the active set are often much more severe than case 1 and 
case 2 changes. 



The effect of a constraint entering or leaving the active set (Cases 1 and 2) can best 
be demonstrated by a simple example from Vanderplaats a id Yoshida (1985). 

Minimize f(x) = 2xi 2 - 2xip + p 2 + 4xi - Ip (2.13) 

subject to: gi=4p + xi>0 (2.14) 

The Lagrangian will be 

L(x,u) = 2xi 2 - 2xip + p 2 + 4xi -4p - ui(4p + xi) (2.15) 

for p = 0, the optimum is f(x*) = 0, xj* = 0, gi =0, and 11 = 4. 

This example will illustrate a constraint leaving the active set (case 2) as p increases. 
The same example can be used to illustrate a constraint entering the active set (case 1) if we 
use a different starting value of p. 

To demonstrate the methods we have talked about, we will calculate the sensitivity 
estimates using four representative methods: the first and second order Kuhn-Tucker 
method, and first and second order extended design space method. We will conclude with 
a comparison of the various methods used to solve the problem. 

To solve for the sensitivity by Kuhn-Tucker equations we use equation (2.6) to 
provide the following system of equations 
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[to 1 ] 

which yield 


dxf 

W 

dui 

LWJ 


= [-4] 


^.= -18 

From equation (2.7) we can 
to the parameter p to be, 

df = <£ + ^L T 2*1 = -4 + 4 (-4) = -20 
dp op 0 x 1 dp 


(2.16) 

(2.17) 

(2.18) 

determine the sensitivity of the objective function with respect 

(2.19) 


The active set will change when the Lagrange multiplier of the constraint goes to 
zero, which can be estimated by equation (2.10) 


A P =^ = ^ =0 - 2222 


( 2 . 20 ) 


(out 

w. 


) 


therefore we ere assured of reasonable results for extrapola ions for which Ap less than 

0 . 2222 . 

For example, a linear approximation by equation (2 8) to estimate the value of the 

new optimum produce 

fnew = f* + Apg| = 0 + Ap(-20) = -20Ap 


A quadratic estimate of the new value of the objecti ve function can be made by 
evaluating the following equation found in Fiacco (1983), McKeown (1980 b), and 
Sobieski and Barthelemy (1983) 


d2f_02L d 2 L 8xi dui 8gL 
d^ 2 _ 3 p 2 + dxidp dp dp dp 

which produces d 2 f/dp 2 = 82. Using the quadratic 


( 2 . 22 ) 

estimat ' for the value of the objective 


function we obtain 

fnew = f* + Ap§ + 0.5ApgpAp - -2«Ap + 41 Ap2 


(2.23) 


The same predictions can be made by Vandetplaa s' extended design space 
algorithm. We begin by formulating the following direction finding problem for decreasing 
values of p, where x 2 represents the parameter p, and x 3 i i an additional vanable to ensure 

that p has the required sign. 

, „ (2.24) 

minimize 4xi - 4x2 * c x 3 
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subject to: xi + 4 x 2 - 0 

-X2 - X3 > 0 
1- (xi2+ X2 2 )> 0 


(2.25) 

(2.26) 
(2.27) 


For c = 1000, the solution is xi = .970142, X 2 - -.212536, X 3 - .242536 which 

yields the following estimates of the sensitivity derivatives 

*20 < 2 - 28 > 

dp 

dxi . (2.29) 

dp -' 4 


For increasing values of p we obtain the following subproblem 
minimize 4xi - 4 x 2 * c x 3 
subject to: xi + 4 x 2 - 0 

X2 - X3 > 0 
1- (xi2 + X22 )> 0 


(2.30) 

(2.31) 

(2.32) 

(2.33) 


When this problem is solved, the resulting sensitivities are sensitive to the value of 
the parameter c. The solution for several values of c are presented below in Table 2.1. 

Table 2.1 The effect of "c" on EDS sensitivity 


variable 

xi 

x 2 

X3 

df/dp 


c =1000 

c=500 

c =100 

c=L 0 

c= 1.0 

O 

II 

O 

b 

-0.398E-2 

0.99999 

0.99999 

-4.016 

-0.79E-2 

0.99996 

0.99996 

-4.0316 

-0.388E-1 

0.99924 

0.99924 

-4.155 

-0.:763 
0.9611 
0.9611 
-5. 50 

-0.624 

0.9611 

0.9611 

-7.196 

-0.707 

0.707 

0.707 

- 8.0 


From this table it is clear that the choice of c will effect the sensitivity derivatives. For 
demonstration purposes c = 10 was chosen, this yielded th e following sensitivity 
derivatives. 


df . 
dp' 
dxi 
dp 


-5.1502 
= -.28756 


(2.34) 

(2.35) 


Vanderplaats and Yoshida (1985) report that the value of c has little effect on the 
EDS algorithm. However Vanderplaats and Cai (1987) re port that after further research the 
value of c will effect the accuracy of the EDS procedure. 


Using Vanderplaats second order extended design space algorithm provides exact 
answers for the sensitivity for this problem. 
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Figure 2.1 illustrates the accuracy of various methods . We can see that when the 
active set changes at Ap = 0.222 the predictions become less accurate. 



the active set 

Figure 2. 1 A plot of various optimal values o f with respect to p 

Figure 2.2 illustrates the location of the optimum vi lue of xi as a function of p, as 
predicted by various algorithms. When the active set changes there is a discontinuity in the 
rate of change of the optimum value of xi with respect to p (i. e. 3xi/3p is discontinuous at 

the point). 
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Figure 2.2 A plot of various optimal values of q with respect to p 


From figures 2.1 and 2.2 it is possible to draw somt conclusions about the relative 
performance of the four different methods that were used tc obtain sensitivity information. 
Using the first order Kuhn - Tucker method we see that the solution follows the inequality 
constraint in both the positive and negative direction. The linear estimate of the new value 
of xi is accurate for small changes in p less than 0.2222. But for values of p greater than 
0.2222, the active set has changed and large errors in the predictions are introduced. This 
is also true for the linear prediction for the value of the objective function. 

The second order Kuhn - Tucker estimate of the value of the objective function is in 
exact agreement in the region where the active set remains the same, as seen in figure 2.1. 
However after the active set changes the predicted value of the objective function is a poor 
predictor of the actual value of the optimum. 

The first order extended design space provides the same results as the first order 
Kuhn-Tucker sensitivity for decreasing values of p. For in creasing values of p we see that 
the search direction changes. This approximation appears o overcome the constraint 
leaving the active set, but it is a poor predictor of the actual value of the optimum for small 
variations in p. For other values of the parameter "c" we v ill obtain similar values for the 
sensitivity derivatives. 

The second order extended design space provides t ie exact values of the locations 
of the optimum value of the objective function. This is because the approximating problem 
that is formulated is the same as the original problem. 

With this simple example we have demonstrated the effect of a constraint leaving 
the active set on the algorithms for estimating parameter sensitivity. We can see from this 
example that, as we might anticipate, using second order estimates can produce more 
accurate extrapolations. In fact, only the second order extended design space algorithm 
provided good results after the constraint left the active set. However its usefulness is 
diminished by the need for second derivatives which can be computationally expensive to 
obtain. 

2.2.4. Example of Case 3 

Recall, that Case 3 is characterized by the adding cf a new constraint to the active 

set and the gradients of the active constraints become linearly dependent. When the 

gradients of the constraints are linearly dependent the Lagrange multipliers will not be 

uniquely determined and the Kuhn-Tucker optimality conditions cannot be uniquely 
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verified. This also results in a discontinuity in the Lagrange multiplier sensitiv d 

When the constraint gradients become linearly depende at for a value of one of the 
parameters i, is assumed dial this is only a temporary condition. If die user is interested m 
die effect of changing the parameter on die optimum then this i iformauon can be obtained 

on either side of the singular point 

This behavior is demonstrated in the following example 
minimize: f = xi 2 + (P - 1) 2 

subject to: gl «3*i + 2P- 10*0 ^ 

g2 = 2 xi + 3 P - 10 > 0 < 2 ' 38) 

When P = 2, the minimum f* = 5 occurs at «• = 2. At this point , both constraints are 
active, and die gradients of the constraints are no, linearly independent. The Lagmnge 

multipliers will be in the family » 

ui,u 2 e {3 ui + 2u 2 = 4,ui >0,u 2 >0} 

A, this point, df*/dp, 3x*/3p and 3u/3p can no, be u, liquely determined. Results 
for these derivatives can be developed if we consider positive and negative changes in P 
separately on either side of this degenerate point which we shall indicate by 3x/3p 
increasing values of p and 3x/3p- for decreasing values of p. 

Figure 2.3 presents the sensitivity plots for this problem. Figure 2.3 (a) and (b) 
represent die firs, order predictions of die new values of the Lagrange multipliers for this 

problem. For this problem the linear P'«“ cd ° ns ^ JT^tional 

multipliers. There is a discontinuity at Ap = 0.0, therefore t e y 

derivatives for these value. Figure 2.3 (c) represents lineai predictions o the new value 

of die objective function. Notice again , dia, "Lot exist for 

prediction and we can not determine df*/dp for Ap - 0. Th..retore P 

this value of p. Figure 2.3 (d) represents the predicted location of x, and we nonce 

same situation as we have for df/dp. 
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Figure 2.3. A comparison of the Sensitivity of a problem 'vdth a Linear Dependence in the 

Constraint normals 


2 3. SUMMARY 

Sensitivity analysis is now routinely used in linear programming (Falk and Fiacco 
1982) and most linear programming algorithms provide modules for the calculation of 
sensitivities. This has not been the case for applications o ' nonlinear programming. The 
most common use of sensitivity derivatives has been in tht area of structural optimization 
and in work done for decomposition methods. Some of the reasons for this may be due to 
a lack of understanding about how to perform sensitivity ; nalysis for nonlinear problems, 
or to a lack of established procedures and supporting software that make the analysis more 
readily available to the average user. The largest contribu or to its lack of use is probably 
the difficulty involved in implementing the current theory and methods. 

An assessment of the methods discussed in Section 2.1 and demonstrated in the 
examples in Section 2.2.4 leads to the following conclusions about the current state of the 
art of parameter sensitivity analysis: 

1 . The Kuhn-Tucker sensitivity equations (2.6) accurately define the desired 

sensitivities assuming no changes in the activ< constraints. To implement these 
equations, however, requires second order information about the Hessian of the 
Lagrangian, and the change in the gradient of he Lagrangian with respect to the 
parameter. Both of which are difficult to obta in reliably for all but a few special 

cases. 
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2 . The Extended Design Space (EDS) method pro' ides sensitivity information 
without the need for the second order information required of the Kuhn-Tucker 
method. However, the sensitivity estimates are effected by a choice of a 
formulating parameter c, and may not give the s ame directions as those obtained 
from the Kuhn-Tucker method. 

3 . Changes in the active constraint set will effect t) le accuracy of any of the 
methods and may limit the region upon which c xtrapolations to the design can 
be relied upon. 
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3. 


New method for Estimating Parameter Sensitivity 


In this chapter a new method for estimating pantmeter sensitivities based on the 
Recursive Quadratic Programming(RQP) method is described. We begin wt a e 
description of the RQP method and the advtuttages it provides for esttmattng — 
Next, we present the RQP based algorithm for estimating para neter sensttivtt) 
the advantages of the RQP method discussed in the previous section. is is o o 

cwuparisonof the new method with existing methods based o 1 l^'llsion 1 

is being produced and the number of function evaluations req rued. Ftnallyj. dtscuss to 
I! rinted of potential problems that may be encounterod wi d, tite new RQP sensmvtty 

method. 


3.1. 


RQP METHODS 


The RQP method has been on the forefront of recent , esearch in optimtzauon 
algorithms and has been emerging as one of the most efftcien t methods available for 
solving small to medium sized, general nonlinear programm ng problems (equatton ^ 
1.6). State of the art RQP methods have been developed by many researchers, sue , 
Powell (1983), Schittkowski (1984), Gill, Murray and Wright (1986) and Bartho omew- 
mggs ( 986, 1987, ,0 name a few. The algonthm has been tested agatns. o*er genera. 
ZL programming algorohms by Schittkowskt (.980), Ecker -^hm, 
,,984). Belegundu and Arora (1985). The results of these , ests have shown the RQP 
method to be one of the most efficient algorithms available 1 or tite solution of nonhnear 

programming problems. 

All RQP methods use the same basic strategy of lint arizing the constrai 
approximating the Hessian of tite Lagrangtan to form a q u,hattc programmtng (QP) 
subproblem The QP subproblem is then solved for the seech dtrecuon, s, and a new 
estimate tite LagrJge multipliers of tite constraints. The QP subproblem has the form 


Minimize l/2sTBs + s T Vf 
subject to Vh T s + h = 0 
Vg T s + g ^ 0 


(3.1) 

(3.2) 

(3.3) 


. „ .hr. Hfsccian of the Laeran 'ian which is normally 

where B is an approximation to the Hessian oi me & , 

Wanted by viable metiic metitod. The lagrange multipliers of tite — for the 
original problem (equations 1. 1-1.6) are estimated by the lagrange multtphers of tite 
constraints in the QP subprob.em (equations 3.1-3.3). Tie search direcuon s ts then used 

to calculate a new estimate of the optimum ^ 


(3.4) 

x it+i = x»t + as 

where a is determined by minimizing a line search penalty function P of the following 
general form, 

P(x,u,v,R) = f(x) + R*n(h,g,u,v) i <35 ’ 

where n represents some combination of the constraints and he Lagrange multipliers. The 

penalty function attempts to assure that both the objective function and the violation of the 
constraints are reduced. As the method converges, the optimal step length a gener y 

approaches 1. 

RQOPT, a typical RQP algorithm, was used in our r< search. A summary of the 
algorithm that is used by RQOPT is presented here, a full description of RQOPT can be 
found in the users manual (Beltracchi and Gabriele 1987 a), or Beltracchi (1985), a e e 
and Beltracchi (1986,1987 b). There were several modifications that were made to RQOPT 
for this work and these will be discussed in section 4. 1 of this report. 


Given x u 

An Approximation to H 
^ and algorithm parameters ^ 

1. Define the Active Set 

I 


2 . Calculate the Gradients and 
update the Hessian Approximation 


3. Solve the QP Subprobler i~| 


1 4. Find the intial step lengti il 


1 5. Conduct the Line Search 1 


j 6. Update Penalty Paramete rs 


Goto Step 1 


Figure 3.1 Flow Chart for RQOPT 

Figure 3.1 shows the basic steps that are used by the RQOPT program. The 
RQOPT algorithm begins with an initial estimate of the lo ation of the optimum and se 
algorithm parameters that have been se, by the user. The nrs, step of the algorithm ts to 
identify the active constraints, it is imponant that the proper constraints are chosen to be in 
the active set as this can effect the rate of convergence of the algorithm and, for our 




purposes, .he approximation of the Hessian of the Lagrangian. Algorithm Peters ™ 
arable to allow the user to control which constraints are considered active dunng 

course of the optimization. 

The second step is to calculate die gradients of die obje. rive function and the 
constraints that are in the active set and then update the approx matron of die Hessian 
Lagrangian. The update of the Hessian is performed using the BFS variable metric u 

with modifications specified by Powell (1977). 

The third step is to solve a quadratic programming subproblem (equatio ns 3 .1-3.3). 
-Hre QP subproblems generated by RQOPT are solved by OPTQP. a spec.nl rmplementaoon 
of the reduced gradient acted. If die subproblem has no tea sible solution^,. tactivese ts 
redefined by dropping constraints from the active set until a ft asible subp 
found. 

The line search for the next point x**' makes up the tourth and frith steps of the 
algorithm An initial step sire for die line search is determined in the fourth step such that 

litis no, in the Jve set are no, excessively violated. "" 

the fifth step, and if a step of o = 1 satisfies the line search a „ena, then that step 

and the line search ended. 

The sixth step updates the penalty parameters used ir the line search, and the 
Lagrange multiplier estimates. We then return to start another tteration. 

Them have been several different variants of die RQP method proposed Some of 
die variants am discussed in Beltracchi (1985). The major differences rn RQP algonthms 

are in the form of the line search objective function (equation 35) and the 

.i / ^ or ill t u at aj-e use d Research continues on these areas 

the QP subproblem (equations 3. 1-3.3) that are usea. 

but no one formulation has yet to prove itself clearly sup ir. 

Some of die penalty functions tha, have been proposed for (3.5) arc a 1, exact 
penalty function (Fletcher 1984, Powell 1987), a 1 2 quadratic loss penalty function 
(Bartholomew-Biggs 1980) or an augmented Lagrangtan (Lhen,Kong an 
1987 ,Barrholomew-Biggs 1985, 1987). The penalty function's parameters am ad usmd 
after each iteration, and how die parameters am updated effects die convergence of die 

method. 

There are two basic philosophies for forming the QP subproblem for ^ 
methods, the inequality constrained (IQP) formulation and equality ’ 

formulation. The most common is the IQP approach which uses a subproblem of the form 
of equations 3.1-3.3. The EQP approach linearizes only t subset e.nequ iy 


constraints and considers these as equality constraints in the subproblem(i.e. equation 3.3 
is considered to be an equality constraint). A discussion of the advantages and 
disadvantages of the IQP and EQP subproblem formulation can be found in (Bartholomew- 
Biggs 1987,1986,1982, Zhou and Mayne 1985, Schittkow ;ki 1983, Murray and Wright 
1982, or Powell 1978). 

Although the method does perform well, it does have some disadvantages. In 
general, the method produces a series of infeasible points v hile approaching the solution 
which may pose a problem for some problem formulations RQP methods are also 
sensitive to variable and objective function scaling and no good scaling algorithms have 
been proposed. Finally, the best penalty function or algorithm for updating the penalty 
parameters for the line search is still a subject of a great deal of research in these methods. 

On the plus side, the following advantages have be on attributed to the method. In 
terms of number of function evaluations, this method appe its to be one of the most 
efficient methods available. This has been demonstrated ii any of the published 
comparison studies in which codes for these methods parti tipants. The method does not 
require a feasible starting point which means there is no sp icial phase 1 search employed as 
in the GRG method or the feasible direction method. Although, as mentioned above, the 
method is sensitive to variable and objective function scali tg, it is not sensitive to 
constraint seating. Finally, the RQP method provides an estimate of the Hessian of the 
Lagrangian, which can be useful for other purposes, and it is very efficient at locating an 
optimum when the starting point is close to the true optimum. Both of these last 
advantages will be exploited in the next section which describes a method for sensitivity 
estimation based on the RQP method. 

3 2. PROPOSED AI .GORTTHM FOR PARAMETER SENSIT IVITY 

In reviewing the current methods for sensitivity analysis in chapter 2, we recall that 
to employ the Kuhn-Tucker sensitivity equations required second order information about 
the Lagrangian. For most engineering problems this type of information is often not 
available in closed form, and estimation techniques would be prone to truncation and 
numerical errors. Therefore, the application of these equations to a broad spectrum of 
engineering applications is limited. 

One proposal mentioned in chapter 2 to circumvent these problems was suggested 
by Armacost and Fiacco (1974) and McKeown (1980 b). Their proposal to estimate the 
sensitivities without estimating the higher order information was given in equations 2.4 and 
2.5, 
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df f(x*'P + Ap) - f(x*, p - Ap) 
dp 2 Ap 


9x* x*(p + Ap) - x*(p - Ap) 

"dp" _ 2 Ap 

These equations represent the use of differencing techniques t o estimate the sensitivines, 
where the values f(x*,p + Ap), x*(p + Ap), etc. are determined by reopommng t e 
problem for the new values of the parameter. For most algor thms, particularly pen ty 
function based methods, the reoptimizations would be a non-mvial task requmng a 
considerable number of function evaluations. However, this is the type of problem where 
the RQP method is considered to be very effective. The goal of the new algorithm ts to 
exploit the strengths of the RQP method to estimate sensitivities by these differencing 

techniques. 


The RQP method possesses two characteristics that we felt can be exploited for 
determining parameter sensitivities: (1) an approximation to die Hessian of the LW" 
is developed, and (2) if this approximation is exact (or close , then the RQP method quickly 
and efficiently solve the reoptimization problem used in the difference equations. 
Essentially, if we can develop g«xi Hessian approximation:., die RQP method is equivalent 
to applying Newton s method to solve the Kuhn-Tucker conditions for the perturbed 
problems which should require only 1 or2iterationsofRQP>. The small number of 

iterations, coupled with the fact that the RQP method should require only a one step hne 
search, should allow the reoptimizations to occur wilhou, ti e need for many function 


evaluations. 


Based on the above arguments, we propose the following procedure to calculate 
parameter sensitivity derivatives (for cases where there are no changes in the active set for 

small variations in the paramters^). 


of the RQP method. . . 

(the * notation is used to denote optimum values) 


Step 1. Perturb the fixed parameter p, to Pi * = Pi c + A Pi where Apt is some small 
perturbation to pi 

Step 2. Perform one complete iteration of the RQP method to find: 
f+ the estimated value of the optimum object ve function 
x+ the estimated value of the optimal of the design variables 


1 We can expect only one or two iter^ions of^RQP if we ;can Should be good. 


directional derivatives. 
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u + the estimated value of the optimum Lagrange multipliers 
g;+ j i Active set 

(as predicted by the RQP method for pi = pi' ) 

Step 3. Perturb the fixed parameter pi to pf = Pi 0 * Ai>i 

Step 4. Perform one complete iteration of the RQP m ethod to find, 
f- the estimated value of the optimum objective unction 
x- the estimated value of the optimal of the design variables 
u- the estimated value of the optimum Lagrange mulhp lers 

gj- j e Active set 

(as predicted by the RQP method for pi = Pf) 

Step 5. Obtain estimates for the sensitivity derivatives from the following central 
difference approximations 

df* f + - f~ (3- 6 ) 

dp"” 2Ap 

dx* x + - x~ (3.7) 

^P"~ 2Ap 

du* u + - u~ (3.8) 

"<5p~“ 2Ap 


Step 6. Estimate the sensitivity of the inactive constraints by 

dgj* _ JLi — I_£l_ j g Active set 
dp 2Ap 


In addition to dte algorithm described above, the fo lowing variants of the basic 
algorithm are also proposed 

1. Forward differencing. For this variant we would omit steps 3 and 4 and then 
use a forward difference approximation (equatir n 3.10 instead of equations 
3.9) to approximate the derivatives 

dq* q+ - q (3.10) 

where q can represent f* , x, u, and the inactive constraints. We may want to 
use this formulation because it requires less function evaluations than the central 
difference approximation. However, the forwa-d difference approximation is 
more susceptible to roundoff and truncation en irs and requires a more accurate 
optimum to yield good sensitivity derivatives. 

2. Forward differences using 2 iterations of the RQP method. This variant is 
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similar to variant 1, but we would perform 2 iterations of RQOPT in step 2. This 
• will yield a more accurate estimate of the optimu m of the perturbed problem. 

When we use this option we can also update the approximation to the Hessian of 
the Lagrangian, or adjust the perturbation Ap to c btain a more accurate estimate 

of the derivatives. 

3. Central differencing using 2 iterations of the RQP method. This variant would 
perform two iterations of RQOPT in steps 2 and 4 of the basic algorithm. As in 
variant 2 we can update the Hessian approximation during each iteration or adjust 
the perturbation Ap to obtain a more accurate estimate of the derivatives. This 
variant is the most computationally expensive ol the proposed variants. 

When there are many parameters that the user need i to obtain sensitivities for then 
the user may want to use variant 2 or variant 3 to calculate the sensitivities for the first few 
parameters. This will allow a more accurate estimate of th s Hessian of the Lagrangian to be 
constructed. After an accurate estimate of the Hessian of tie Lagrangian is built, the user 
should switch to either the baseline or variant 1 to obtain the sensitivities of the remaining 
parameters. The Kuhn-Tucker sensitivity equations may i Iso be used with the Hessian 
approximation, after a good estimate of the Hessian of the Lagrangian is built. However 
the Kuhn-Tucker sensitivity equations also require 3V x L/c p be calculated and this term 
may be subject to numerical noise because V X L = 0. 

3 3. COMP A R ISON TO EXISTING METHODS 

This section provides a derivation that indicates the performance that is expected 
from the new sensitivity algorithm. This section also presents a comparison between the 
RQP based method and two existing methods described ir chapter 2 based on the number 
of function evaluations required to estimate the sensitivities. 

3.3.1. Demonstration of Equivalence of New Me thod to Kuhn-Tucker Metho d 

This section will show that the finite difference approximations obtained by the 
proposed method are in fact equivalent to the sensitivities obtained by solving a modified 
set of Kuhn-Tucker sensitivity equations. The modification of the Kuhn-Tucker sensitivity 
equations involves replacing the Hessian of the Lagrangi; n with the approximation B, 
obtained from the RQP method. 

The following assumptions are made for this derivation; no equality constraints are 
present, the base optimal point is stable 3 , and the gradient are continuous. The derivation 

3 A stable point is defined as a point where the acitve set does not i hange for small variations in the 
parameters 
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in the presence of equality constraints does not change too much but the equality constraints 
were left out to simplify the notation. If the base point is not stable then this derivation can 
be used to find directional derivatives; this will be discussed at the end of this section. If 
the gradients are not continuous then we cannot even be assured of an optimum point since 
the assumption of continuity is also made for the derivatior by the Kuhn-Tucker method. 


We begin by restating the Kuhn-Tucker Sensitivity equations 


r 0 

dx 


av x L _ 
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We strive in this derivation to show that the proposed method is equivalent to estimating the 
sensitivities using modified version of equation (3.1 1) that replaces V^L with B obtained 

from the RQP method. If this is the case, then we can anticipate the kind of accuracy to 
expect and where the possible sources of error will result. 


If we examine the equations (3.6-3. 8), used by the proposed RQP sensitivity 
method we see that these provide finite difference approximations to the sensitivity 
derivatives of the objective function, design variables, and Lagrange multipliers with 
respect to pj. The derivatives are defined by the following 


df* _ lim 
dp - Ap->0 


/ T*(x*+Ax,p°+Ap) - f*(x*,pOn 
, A P > 


( 3 . 12 ) 


dx* _ lim 7x*(p°+Ap) - x*(p°) 
dp Ap->0 ^ Ap 

A*=^-Ap 


3u* lim 
Ap— *0 


^u*(pO+Ap) - u*(p°)^ 
< A P J 


where p° represents our base point. 


(3.13) 

(3.14) 

(3.15) 


The RQP subproblem for the simplified case wher i the active constraints remain 
active and there are no equality constraints can be written is 

min 1/2 s T Bs + s T V x f (3.16) 

subject to s T V x gj + gj = 0 j € Active Set (3.17) 

where B is the approximation to the Hessian of the Lagrangian and the inequality 
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constraints gj are considered as equality constraints. 

If we assume 4 a step length of ex = 1 is used in the 1 me search (equation 3.4) we 
can rewrite equation 3.4 in terms of x' - x as 

s = X' - x < 318 > 

where x’ is the new estimate of x*. Substituting equation 5.18 into equations 3.16 and 
3.17 we obtain the following subproblem which is minimized with (x - x) as the design 
variables 


min 1/2 (x’-x) T B(x'-x) + (x'-x) T V x f(x,p + ApO 


(3.19) 


subject to (x'-x) T V x gj(x,p + A Pi ) + gj(x,p + ApO =0 j e Active Set (3.20) 
We can now state the optimality conditions for the subproblem represented in equations 
(3.19-20) as 

B(x’-x) + V x f(x,p + Api) - u'V x gj(x,p + Api) = ) (3-21) 

(x’-x) T V x gj(x,p + A Pi ) + gj(x,p + A Pi ) = 0 j € Active Set (3.22) 


Here u' represents the estimated value of the Lagrange mu tipliers at the new optimum. 


Now we substitute into equation (3.21) the following definitions of zero 


V x L(x,p°) = 0 = V x f(x,p°) - uV x g(x,p°) 

uV x g(x,p°+ Ap) - uV x g(x,p°+ Ap) = 0 
This will yield 


(3.23) 

(3.24) 


B(x'-x) + V x f(x,p° + Ap) - u'V xgj (x,pO+ Ap) . ;V x f(x,p°), uV x g(x,pO))+ 
uV x g(x,p°+ Ap) - uV x g(x,p°+ Ap) = 0 (3-2: 

Rearranging we obtain 


B(x'-x) - u'V x gj(x,p°+ Ap) + uV x g(x,p°+ Ap) + 

(V x f(x,p° + Ap) - uV x g(x,p° + Ap)) - (V x f(x,p°) - uV x g(x,p°)) = 0 (3.26) 

Rearranging further and writing in terms of the Lagrangia n function we obtain 

B(x'-x) - (u' - u)V x gj(x,p°+ Ap) + V x L(x,u,p° t- Ap) - V x L(x,u,p°) = 0 

(3.27) 


Now we will divide equation (3.27) by Ap and take the limit as Ap goes to zero to 


4 A common assumption for RQP methods 
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obtain 

to RtV-vi ■ III' - u)V,g j lx.D°+ M * v.u»-u.p°.±.>p> - v » L( ’ t - u ' ^ = 0 
4p - 0 iP (3.28) 

Using .he additive and multiplicative properties of the Limit fur cion we obtain 
_lim /x'-xn lim 11111 n (V xgj (x,p°+ Ap)) + 

B ap-*o(a P ) Ap "°l A P J P 

li m |^ VxL(x,u,p°+ Ap) - VxL(x,u>P^r |_r> (3.29) 

Ap— >0 1 Ap J 

Now we can use the definition of a derivadve of some function h with tespec, to some variable p 


3h _ lim h(p+Ap)-h(p) 
<5p Ap— >0 Ap 


(3.30) 


Applying the definition of 3x/dp, 8u/dp to (3.29) we obtain 

-v -v ,. 9V x L(x,u,p°) _ n (3 31) 

B| - | a"o V ^ (X ’ P ° +AP) + 

If we use the standard assumption that the functions are twice continuously differentiable 
we can state 

(3 32) 

^0 Vxgj(x,p0 + Ap) = Vx8j(X,p0) 

And notTsubstituting equation (3.32) into equation (3.31) we obtain 

B& - V xgj («,p0)| + 5Zi|^=0 < 3 - 33 > 

The equation above (equation 3.33) represents the firs, par, of 

equations with tire approximation B instead of the Hesstan , ,f the Lagrang . 

The next step in this derivation is to examine equation (3.22) in terms of P° ♦ Ap 

we obtain 

(x’-x) T V x gj(x,p° + Ap) + gj(x,p° + Ap) = 0 j e Active Set 
Now we can subtract gj (x,p) = 0 from equation (3.34) to obtain 

(x’-x) T V x gj(x,p° + Ap) + gj(x,p° + Ap) - gj(x,p°)= 0 j e Active Set (3 

If we divide equation (3.35) by Ap and take the limit as A, goes to zero we can write 
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lim 7(x'-x) T V x g i (x,p 0 + Ap) t g j(x,p° + Ap) - gj(x,p°)^ _ Q ^ ^ 

Ap->0 Ap Ap J 


Using the additive and multiplicative properties of the limit function we obtain 

. n . - lim /x'-xn , Urn (z (x,p° + Ap) - gj(x,p°)A 

, 0 + a p> * a^ 0 (i7) + ^{- — J- 0 


lim 

Ap- 


(3.37) 


Again using the definition of a partial derivative of (equati in 3.30) and we obtain from 
equation (3.37) 

dx d 


!£* V xgj(x,p° + Ap)»f 
Using the results in equation (3.32) we obtain 

Vxgj(x,p°)*|^ +^=0 

Which represents the second part of the Kuhn-Tucker sensitivity equations. 


(3.38) 


(3.39) 


Now equations (3.33) and (3.39) can be assembled into matrix form to yield 




dx~ 


"av x L" 

B 

-V,gT 

-V x g' 

0 

du 

L¥-J 

+ 

Li J 


= 0 


(3.40) 


Equation 3.40 is the same as equation 3.19 with the excep tion that equation 3.40 uses the 
approximation, B, of the Hessian of the Lagrangian in pi; ce of, V^L, the true Hessian of 

the Lagrangian. Referring to (3.40) as the modified Kuhr -Tucker equations, we see that 
the proposed method is principally a difference approximation to the modified Kuhn- 
Tucker equations. This implies that if B is a good approximation of the Hessian of the 
Lagrangian, and a proper choice can be made for the difference parameter that minimizes 
truncation and roundoff errors, then we can produce sensitivity derivatives without the 
need to obtain or estimate the second derivatives required of the Kuhn-Tucker method. 


Several examples were tested to see if the sensitive ty derivatives estimated by the 
RQP method with one iteration converged to the value of sensitivity derivatives estimated 
by the Kuhn-Tucker sensitivity with the approximate Hessian. From these examples we 
observed that the sensitivity derivatives delivered by the rew RQP algorithm are close to 
the derivatives approximated by the Kuhn-Tucker method with the Hessian approximation 
One of these examples is presented here to show this agreement. 
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Test problem 2 (which is described in ihe appendix) is u ted to demonsbaie Ita 

0 m U)' sol v" die problem in 
one iteration and yielded .he following approximation ,o the He ssian mainx* 


H approx — ; 


3 2 2 1 
2 3 2 
.2 2 3 J 


— OFF*'"' 9 9 3 J 

„ we use this Hessian approximation to solve for the sensitivi, , of paramemr 1 by et^anon 
(3.40) we will obtain the following system of equations 


■3 2 2 -1 1 
2 3 2 -1 
2 2 3 -1 
.1110 


3xTl 

dpi 

3x2 
c5pT 
3x3 

3pT 

i 9ui i 

L-SpT-l 


+ 


r-12 
5 
0 
1 


= 0 


(3.41) 


the solution of these equations yields 

_ (9 33333,-7.66666,-2.66666) 

dpi 

9u 


Bpi 


= (-4.66666) 


(3.42) 


(3.43) 


The RQP based sensitivity algorithm calculated the following sensitivity derivative 
approximations. 

(3.44) 


3 * . = (9 33333,-7.66666,-2.66666) 
dpi 

3ui 


(3.45) 

= (-4.66666) 

The above derivatives were calculated using the RQSEN program (described in section 4 

rr^ru, 

differencing, equations 3.7 ,3.8) ana one ucia 

,f tire base point, p», is unstable (degenerate) we o use a similar derivation to 
If the base point, P , f , for Acting the sensravtues of the 

calculate dimaiona. derivauves, whtc wdl be u^ful o , wU1 * 

design variables and Lagrange mulupliers. The use of dm cuonai 

discussed in section 6. 

SThe Hessian approximation for problem 2 is 

appendriof thisVt). ™s is t because the comparison to the Kuhn-Tucker 

w^coold clearly indicate the perfonnmee “gp/ 
sensitivity method with the approximate Hessian trom W 
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3.3.2. Performance Comparison with Other Methc ds 

This section compares the RQP based method to tw o of the methods discussed in 
chapter 2; the Kuhn-Tucker method, and the extended design space (EDS) method. The 
comparison is based on Table 3.1 which examines the number of function evaluations 
required by each method to calculate parameter sensitivitie 5 df*/dp, dx*/dp and du*/dp 
(assuming that when the optimum is found that the Kuhn-3 ticker conditions have been 
checked, this means that V X L and the Lagrange multiplier are known before the 
sensitivity analysis is performed). It is assumed that the ob jective function and constraints 
are interrelated 6 . It is also assumed that problem linearity or problem form are not 
exploited in calculating parameter sensitivities. 

The first row of Table 3.1 represents the methods used in this comparison. The 
second row represents the number of function evaluations required to calculate the 
sensitivity derivatives for the first parameter. Subsequent parameters may require fewer 
evaluations for some methods. 

The first column of Table 3.1 represents the numtx r of variables present in the 
problem. The second column represents the amount of wc >rk required to solve for the 
sensitivities using the Kuhn-Tucker sensitivity equations. The third and fourth columns 
represent the number of function evaluations required by the EDS algorithm. Column 3 
represents the first order method and column 4 the second order method. The fourth 
column, RQP 1, indicates that forward difference approximations were used to calculate the 
gradients. The fifth column RQP 2, also uses forward difference approximations but 2 
iterations of RQOPT are allowed during the reoptimizatioi . The fifth column RQP 3 
represents the amount of work required for the base line al gorithm using central difference 
approximations. The sixth column RQP 4 represents usin g central difference 
approximations with 2 iterations of RQOPT. 

If the objective function sensitivity is calculated bj equation 2.12 
(df*/dp = df/dp - u dg/dp) then assuming that objective i nd constraint information can 
both be obtain in one call, only one extra function evaluation is required to determine df/dp 
and dg/dp. However, if one wants the design variable anc Lagrange multiplier sensitivity 
then some other equations must be used. 


6 The value of the objective function and all of the constraints are ca culated by one subroutine 
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n Kuhn-Tucker EDS (1 st ) EDS ( 2 nd) 
n2 3n , . (n+l)2 3(n+l) 

T + T +l 1 2 + 2 


RQP 1 

RQP 2 

RQP 3 

RQP 4 

n+1 

2 i+2 

2n+l 

4n+2 


1 

4 

1 

5 

2 

6 

1 

9 

3 

10 

1 

14 

4 

15 

1 

20 

5 

21 

1 

27 

10 

66 

1 

77 

15 

136 

1 

152 

20 

231 

1 

252 

40 

861 

1 

902 


2 

4 

3 

6 

3 

6 

5 

10 

4 

8 

7 

14 

5 

10 

9 

18 

6 

12 

11 

22 

11 

22 

21 

42 

16 

32 

31 

62 

21 

42 

41 

82 

41 

82 

81 

162 


Note: RQP 1 uses forward difference approximations and one item don to solve the perturbed problem 

RQP 2 uses forward difference approximations and two iterations to solve the perturbed problem 
RQP 3 uses central difference approximations and one iteration to solve the perturbed problem 
RQP 4 uses central difference approximations and two iterat tons to solve the perturbed problem 


Table 3.1 Comparison of Various Algorithms for tse in sensitivity analysis 


The following observations can be drawn from this table. 

1. For the Kuhn-Tucker sensitivity equations, most of the work in finding the 
parameter sensitivity is involved in the calculation (by finite differences) of the 
Hessian of the Lagrangian. However, after the first parameter sensitivity is 
determined the cost of evaluating successive ser sitivity derivatives is reduced to 
(n+1) extra function evaluations. 

2. For the first order EDS algorithm, the work required to calculate the parameter 
sensitivity does not increase with problem size. However, this algorithm will 
not deliver du/dp and this algorithm may not be able to find the correct value for 
3x/dp. This will mean that df*/dp will also be inaccurate with this method. If 
the problem is fully constrained the accuracy of dx/9p is better but the method 
may still provide inaccurate derivatives. 

3. For the second order EDS algorithm most of th ; work is in the calculation of the 
Hessian of the objective function and the Hessi in of the constraints. The work 
involved for calculation of successive parameter sensitivities only requires 
approximately n+2 extra function evaluations, rhis algorithm requires the 
solution of a quadratic approximating problem or every new value of the 
parameter supplied by the user. 

4. For the RQP 1 algorithm (forward differencing) is the most efficient of the RQP 
methods proposed and seems to be much more efficient than the Kuhn-Tucker 
algorithm. The work required to calculate successive parameter derivatives is 
constant (n+1 function evaluations). This algorithm will perform well when B is 

a good approximation and the perturbation Ap s properly chosen. 

5 The number of function evaluations for the RQP 2 algorithm (forward 

differencing and 2 iterations of the RQOPT algorithm) grows linearly. The work 
required to calculate successive parameter derivatives is constant (2n+2 function 
evaluations). The work for calculating success ve parameter sensitivities may be 
reduced because the Hessian approximation wi l improve after each parameter 
sensitivity derivative is approximated, which v. ill eventually reduce the amount 
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of work required to solve the perturbed problem. 

6. For ihe RQP 3 algorithm (central (2^ 

and the work for calculating successi P g . Qf ^ sens itivity derivatives 

<•*»*» ° f *• " s as follows 

(j2f f+ - 2f* + f' (3.46) 

dp 2 Ap 2 

-This approximarion of the second ^ 

a^mage^f^sii^'wnn^dif^ences occurs when rite active set changes and 
directional derivatives can be approximated. 

7. The RQP algorithm with central ddferen^g and . ^^or^r^ui^t^calculate 
most expensive of the proposed The 

SSSTSfp iterations Seas 

approximation is improved. 

The above discussion dealt with me number of requited htncrion evaluations to 
calculate the parameter sensitivities. We did not account fot any of the otiter overhead su 
as solving the QP subproblem fot the RQP method or solving a quadratic approxtma g 
problem for the second order EDS algorithm. 

The overhead associated with using the Kuhn-Tuckcr sensitivity equations is 
relatively small after the fust parameter sensitivity is calculated, tins ts because = > a 
factorization (i.e. LU) is used to solve the Kuhn-Tucker sers.tiv.ty equations then the 
"of overhead becontes o(n, Hops. The overhead fo, solving the t RQP subpro lems 
Tui also be realitively small if a good implemenuttion of the RQP method ts used (re a 

fust order EDS metftod will also be relatively small. However the overhead for the second 
order EDS method could be large depending on the problem. 

In summary, the RQP based methods are competi. ve with the existing methods 
AU variants of the RQP based methcxl requue approximate ly the same number o unction 
evaluations for small problems (n<5), but considerably 1. ss for larger problems (t»5). 


POTF.NTIAI f PROBLEMS 


14 _ 

One of the main issues that needs to be investigate concerns the Hessian 

approximation: will the approximation converge in practice as predicted by e eory. 
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convergence has not taken place then we need to investigate how to improve the Hessian 
approximation. Some modifications that can be made to ol >tain a more accurate Hessian 
approximation are discussed in Chapters 4 and 5. 

As with the estimation of any gradient by finite dif erences, the perturbation step 
size Ap and the nonlinearity of the problem will effect accuracy of the derivative 
approximation. Rules from Gill, Murray and Wright (1983) or Adelman, Haftka, and Iott 
(1986) can be investigated as a means to select the step siz ; Ap. An automated selection 
proceedure for Ap should be investigated after the initial RQP sensitivity algorithm is 
tested. 


When using the forward difference option the choi ce of Ap is even more critical. If 
Ap is too small and the optimum of the problem is not known exactly then when the 
perturbed problem is solved we may only be seeing a bettor estimate of x* being found 
rather than an estimate of the solution of the perturbed pro :>lem. This will cause the 
derivative approximations to be inaccurate. If Ap is too la 'ge then we may only be 
obtaining trend information for the problem. 

All optimization programs incorporate some kind of convergence criteria that is 
based on the relative change in the design variables. This stopping criteria will effect the 
calculation of the sensitivity derivatives for all available methods, because there is a 
common assumption that the base point is a true optimum The central difference 
approximation may be less sensitive to inexact solutions because the solution of the 
perturbed problems will be of a similar degree of accuracy. 

When solving the quadratic programming subproblem some type of convergence 
criteria is normally used. How small this tolerance is will effect how much work is needed 
to solve the subproblem (Nash 1985). During the early stages of the optimization it is not 
advisable to locate the exact solution of the QP subproble n as this may be too expensive. 
However once the program is in the region of a minimum the solution of the subproblem 
needs to be accurate. Therefore, we expect to use a tight convergence criteria for our QP 
solver during our reoptimizations. 


3.5. SUMMARY 

We have proposed a method and some variants b< sed on the RQP method for 
estimating parameter sensitivities which provides sensitiv ty estimates nearly equivalent to 
the Kuhn-Tucker method. The method avoids the need for calculating second derivatives 
and its efficiency is competitive with current methods. The accuracy of the method 
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depends on two major pieces of informa.ion, the quality of .he Hessian approximation 
provided by the RQP meihod, and the step size of the difference parameter us in 
difference formula. Bod, these aspects of the method will be iiscussed m the followtng 


chapters. 
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4 . Implementation 


This chapter will discuss the implementation of the ne w parameter sensitivity 
method described in chapter 3. The program used as the basi > for testing the new method 
was the RQOPT program which is an implementation of an at tive set RQP method 
(Beltracchi and Gabriele, 1987). The discussion begins with a discussion of the 
modifications made to RQOPT to perform the necessary calculations, and ends with a 
description of the software system developed to calculate pan .meter sensitivities. 

4.1. MODIFICATIONS TO RQOPT 

Most of the modifications to RQOPT were concentrat id in one of the major areas of 
concern for the new sensitivity algorithm, the Hessian approximation. These modifications 
are discussed in subsections 4.1.1 and 4.1.2. The line search of RQOPT was also 
modified to yield a smoother convergence to the problem solution and this is discussed in 
subsection 4.1.3. The final modification discussed in subsection 4.1.4, provided the 
option of using a different variable metric update to yield a more accurate Hessian 
approximation 

4 11 Implementation of a Factorized BFS Variable Metric Up date 

Variable metric updates have been successfully used for the past 20 years for 
unconstrained optimization and have been used successfully for approximately the past 10 
years for constrained optimization. Variable metric updates t ttempt to build an 
approximation to the Hessian matrix using only first order information, and solve for the 
search direction from the following equation 

s = B- 1 Vf C 4 - 1 ) 

where B represents the approximation to the Hessian, Vf the gradient of the objective 
function, and s the search direction of the design variables. V; ariable metric updates have 
been provided in the literature for approximating either the inverse of the Hessian or the 
Hessian itself. 

Variable metric updates all have the same basic form. They begin with an 
approximation to the Hessian matrix, and then update the approximation by some rank one 
or rank two correction. The form of the update is normally 

Knew = Bold + vv T + ww T ( 4 - 2 ) 

where v and w are calculated as some product of the old He ssian approximation, the last 

search direction, and the change in the gradient of the objective function. 

Several different forms of equation 4.2 have been p oposed. The most popular 
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variable metric update has been the BFS (also known as th a BFGS) which was proposed in 
1970 simultaneously by Broyden, Fletcher, Shanno and Goldfarb. The BFS update has 
been shown to be the best general purpose variable metric i ipdate. 

One of the problems associated with the BFS variat le metric update is that it is 
effected by the problem scaling. Shanno and Phua (1978) have proposed a self scaling 
version of the BFS update. Its use in a RQP algorithm wa: investigated by Van der Hoek 
(1980). He found the self scaling variant with the second Oren-Spedicato (Oren 1974) 
switch seemed to perform the best with the particular RQP algorithm that he was using. 

In the mid 1970's several authors proposed updating the LDLT factors of the 
Hessian approximation with a procedure that could be usee to stabilize the BFS update in 
terms of the numerical noise encountered in the calculation of the update. With the LDLT 
update we can be assured that the Hessian approximation remains positive definite, this will 
assure that the search directions that are generated from (4 1) are downhill. Additionally, 
finding the search direction from equation 4.1 becomes a s mple matrix calculation when 
using the LDL t update 

When variable metric updates are used for RQP methods it is normally preferred 
that the approximation of the Hessian of the Lagrangian be updated instead of its inverse. 
This is because solution of the QP subproblem requires th ; Hessian approximation. The 
BFS variable metric update is used by most of the success ill implementations of the RQP 
method. 


The BFS update that was used in RQOPT is defined as 

?T 


®new — Bold 
where z and w are defined as 


(zBoid )(zB 0 id) T w w 
‘ zTBoidZ 


w T y 


z — x new " x old 


y = V x L(x n ew,Vnew>Une\v) - V x L(x 0 ld,Vnew»U n ew 1 


0H 


l 

0.8 z T Bz 


[z T Bz - z T y 


if z T y > 0.2 z T B z 
otherwise 


(4.3) 


(4.4) 

(4.5) 

(4.6) 


w = © y + (1 - @)Bz v*-') 

Where the 0 term in equation 4.6 and 4.7 was defined b\ Powell (1977) to help maintain 
positive definiteness of the Hessian approximation, under normal operation 0 is equal to 
The Hessian approximation is guaranteed to be pos tive definite if z T w is greater than 
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zero. The Hessian approximation is no. updated by RQOPT if z T w is less than zero. 

For this study, the LDLT update for the BFS variable netric (defined in equation 
4.3) as described by Gill and Murray (1978) was implement, d (where z and w were 
calculated by equations 4.5 and 4.7). This update uses sever; d matrix transformations to 
achieve a stable update. The actual update of the Hessian approximation is performed with 
a procedure described by Fletcher and PoweU (1974) and ex, ended by Gtll. Murray, and 

Saunders (1975). 

In addition to the stability of this update relative to m unerical noise, as discussed 
above the LDU update provides a convenient means for establishing a reset criteria for the 
Hessian approximation. The need for a reset of the Hessian approximation is discussed in 

the following section. 

4 1 2 Condition Number Reset 

Occasionally, due to numerical noise or a highly nor Unear problem, the Hessian 
approximation may become singular or indeftnite. When ih is happens we can no longer be 
certain that the resulting search directions will satisfy the descent property that ts assumed 
by the RQP The only means to recover from this situation is to reset the approximation to 
some known positive definite matrix, which is generally the identity matrix. Early version 
of the BFS update were reset every n + l iterations but this it a conservative approach that 
will sometimes erase good information and slow the convergence of the algonthm. The 
current thinking is to use a less conservative reset criteria tl, at is based on a condition 
number estimate of die matrix with the hope that useful information built up in prevtous 
iterations is used for more iterations and should result in better convergence. 

The original version of RQOPT reset the Hessian approximation every time the 
active set changed or every n + l iteration. A change in the active set results in a different 
QP subproblem to be solved and it was felt that the Hessia t approximation would no 
longer be valid. Using this conservative reset criteria wood prove unacceptable if we were 
using RQOPT to perform sensitivity analysis. With this reset criteria, we nsk resetting the 
Hessian approximation just before the optimum is reached and would be left wtth only a 
few iterations of the method upon which to build an approximation. Thus we may have a 
very poor Hessian approximation when it comes time to perform the sensitivity an ysts. 

The reset criteria adopted has been used successfully by several other algorithms 
(Powell 1985, Schittkowsld 1983, Arora and Tseng 1987 ). The new reset criteria resets 
die Hessian approximation when the estimate of the cond tion number exceeds a fixed limit. 
This estimate can be found by computing 
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(4.8) 


cond(H)est = ^S' nlT 

where dmin and dmax are ^ smallest and largest values of the D natrix in Ihe LDL 

factorization. 

■ „• criteria has led to a more stable update yielding faster convergence 

for the RQOPT pre^rnand more accurate estimates of the Hes rian of the Lagrangtan. 


Kyuri \n ug^m***** — — 

4 ! 3 r- a ln..laQ. ril"- 1 ^grange h'lnltiplirr F.snmaE S 


The Lagrange multiplier estimates are an reused as hearts tiT the 

appmximation. The value of * e Ugra ,,gian function, 

variable metric update to approxrmate the Hesstan 

The original version of RQOPT calculated .he UPW 
Lagrange multipliers of the constraints ^ a step of a = 1 is used in 

to the true Ugtange multipliers as the problem converges. 

in the first few iterations of RQOPT. At the 
A problem can arise, however, in the ;st]mates produ ced by the QP 

beginning of a search it is possible that a 8™* ^ value of the Lagrange 

subproblem will be several orders of magnitude g ^ ^ value „f die 

multiplier. If the line search then mate a . ^ ^ ^ appK)ximation in such a 

Lagrange multiplier estimate may ^ < ^V m witfl the large Lagrange 

corrected. 

RQOPT was modified to use die following linear interpolation to update die value 
of die Lagrange multiplier estimates after die line search is completed ^ 

, n r of :°“ is'usedin "t search Ration 3.4, then formula 4.9 
updates Lagrange multipher estimates to be the estimate s : riehvered by the QP 

subproblem. This update was also used by Schittkows . 

updates. 


40 



4.1.4. SRI update 


The SRI update is a variable metric update that do is not require exact line searches 
for quadratic convergence, where as the BFS update requires exact line searches for 
quadratic convergence. Because the RQP method seldom performs exact line searches, it 
was felt the SRI update may be able to obtain a better approximation of the Hessian of the 
Lagrangian. 

A table describing the differences between the BR and SRI update is presented 

below 

Advantages Disadvantages 

update 

BFS Self Correcting Requires exact line searches 

Stable (maintains positive definiteness) 

Has a good performance history 

SRI Does not require exact line searches upxiate may be undefined and it is not 

guaranteed to maintain positive definiteness of 
the Hessian Approximation. There is not alot of 
literature on the performance of this update. 

Table 4.1 A comparison of the BFS and SRI variable metric updates 

The stability of the BFS variable metric upxiate has led to its use in almost all RQP 
implementations. However Cha and Mayne (1987) report Jiat they have tested the SRI 
update and found exact convergence of Hessian approximations for quadratic functions. 

Although the SRI upxiate lacks the stability of the BFS upxiate, we were interested in 
comparing the performance of the 2 upxiates in terms of the Hessian convergence. If the 
SRI upxiate delivers better Hessian approximations than the BFS update then we will have 
to further investigate methods to stabilize the SRI update. 

The SRI update is defined as follows 

» d . . i. ( R oidy - z)(B 0 idy - z) T ..... 

new ■ old yt(Boidy • Z) (4 ' l0) 

where y and z are obtained from equation 4.4 and 4.5. This update is undefined when the 
denominator is equal to zero. The SRI update may be undefined even for positive definite 
quadratic problems. This problem was addressed by BraytDn and Cullem (1979), Cullem 
and Brayton (1979). 

The symmetric rank one (SRI) update was implemented in both a factored (LDL T ) 
and unfactored form. In our implementation if the absolute value of the denominator (in 
equation 4.10) is less than some small number we use the HFS update which is described 
in section 4.1.1. 
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Even though the SRI update may be undefined, it has he very nice property of not 
requiring exact tine searches. This is imporian, because in the RQP method we do not 
perform exact line searches, and the BFS variable metric method assumes exact 1 
Inches Powell (.986) clearly demonstiates die detrimental effect of inexact hue searches 
on dm BFS method. The performance of die SRI update for giving quadrauc problems ts 
such that after n updates (providing that all updates arc define 1) die Hesstan approx.ma 
wd. have converged to the hue Hessian. Thus we may obmin a beder convergence of the 
approximation of the Hessian of the Lagrangian if we are abl. = to use the SRI update. 

Some preliminary results were obtained comparing the BFS and SRI variable 

metric updates and these are discussed in section 5.3. 

iMATICALLV CALCULATE 



In this section we provide a brief overview of the sol tware system created for 
studying parameter sensitivities. The software system is made up of three major pieces, a 
problem preparation package RQCRE, the RQP algorithm using dm ” 0<h " 
described in the preceding section, RQOPT, and an interact, re program RQSEN> 
as a post processor/sensitivity analysis module for die RQOPT program me R(^PT 
’ was an existing program and has *en documented previously (Be.tmccht and 
Gabriele. 1986). The RQCRE and RQSEN programs were created I fortius ;« Ml 
be briefly described in the following paragraphs. A more detailed discuss 
systems is provided in the appendix. 

42A The ROCRF Support System 

The RQCRE program is set up to be used as an interactive tool for use with the 
RQSEN system. The purpose of the RQCRE program is to remove ^chance oferro.m 
the pmblem formulation. Die RQSEN program requires approximately: » amty 
dimensioned which are automatically dimensions by 

automatically writes the calling program and dam files required by die RQSEN syste . 

The RQCRE program requires the user to provide basic information about dm 
pmblem such as the number of variables, number of equality constraints, number of 
inequality constraint, and number of parameters that while studied. 

The RQCRE system then produces a main calling program, a shell of the func ' 
subprogram 1 used to define the objective function and the constraints, and a dam file used 
for input into the RQSEN system (sample output is provided in the appen tx). e 
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RQCRE program also sets up the default values for the algorithm papers used by 
RQOPT. 

I" Ths pn *EN nroeram 

The RQSEN system was set up as a pie and R ^ ram for 

program. Tta RQSEN system was set up to be an tntemcttve user frtendly p gram 

performing the following basic functions; 

1. The system can be used to solve optimal design problems 

2 The system can be used to calculate parameter sensib vities 

3. The system can be used to conduct studies of large variations in problem 

parameters ^ ^ ^ ^ IO CTeatc sens.ttvtty plots t f tita. can be used to 
perform trade off studies. 

A sample session with the RQSEN system illustrating these o, ttons is presented 
appendix. 

The RQSEN system requires a calling program and a uncnon^^gram ^ 

(defining tire objective function and the algo nthm 

RQSEN system a,s„ requires titeu^« ^ ^ , tesign param eters. The user 

^“"oSPN^program to study the sensitivities of only cenam parame K rs. 

The RQSEN program first produces optimum ^ ^j^ves, 

optimized the RQSEN system can be used to ' vamtions in the 

which can tiren be used to srndy the effect on tire °P<™ m « program can 

parameters. The RQSEN system is also up s ^ , he paiametm can be 

used to create plots of the optimum sens.nvtt.es for large 
studied. A typical plot is presented in figure 

an aid in creating tire calling program and function subprogram. 
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+ f(X*) 


P,„s similar ,o this one can also be generated for the design variables, Lagrange 


multipliers and v ^j* es 

of chapter 5. Plots similar to figure 4.1 are presented m the . ppendut 


test set 
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5. 


Numerical Experiments 


This chapter describes the numerical experiments t hat have been conducted to date 
on the new sensitivity method. We begin by discussing the initial test set used and any 
special features of the selected problems. Next, we discuss testing that has been performed 
comparing the accuracy of the known Hessian to the approximations obtained, which 
includes comparisons of the BFS and SRI updates. In the third section, the accuracy of 
the sensitivity derivatives obtained with the new sensitivity algorithm is assessed against 
the known results. This section also compares the effect of choosing a central or forward 
difference formula and the effect of the step size Ap. The final section presents some 
conclusions drawn from this initial testing. 

5JL INITIAL TEST SET 


A two phase testing program has been formulated for studying the effectiveness of 
the new method for estimating parameter sensitivity. The first phase was to develop a set 
of test problems for which the parameter sensitivities coul d be exactly determined using the 
Kuhn-Tucker equations. This required that any second Older information needed could be 
determined analytically. Choosing problems of this type would allow a direct comparison 
of the sensitivity results produced by the new method with the exact sensitivities and also 
allow the comparison of the BFS and SRI Hessian approximations. From this study we 
hope to develop some insight into several questions concerning the algorithm such as: 
proper choices for algorithm parameters (i.e. the proper size of Ap), what is the most 
reliable Hessian approximation, how close does the Hessian approximation have to be to 
achieve good results, does updating the Hessian approxin lation during the sensitivity 
analysis significantly improve the estimate, and which of :he variants (forward/central 
difference approximations with one or two iterations of RQOPT) described in chapter 3 
provides the most consistent results. 

The second phase of the testing would consist of testing the algorithm against a set 
of engineering problems where second order information would not be available. Here the 
results obtained from the sensitivity algorithm would be compared to actual reoptimization 
results to assess its accuracy. In the time allotted for this study, only the first phase of 
testing has been completed and is reported on here. 

The problems making up the initial test set are presented in the appendix of this 
report. We have experimented so far with 4 test problems that have a total of 12 
parameters. The problems possess both linear and nonlinear behavior. We expect to 
expand this test set in the near future. Plots of the optimum sensitivity for selected 
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problems and parameters are also presented in the appencix. 


5.2. CONVERGENCE OF THE HESSIAN APPROX 1MATION 

The derivation given in section 3.2.1 showed that equivalence of the new method 
with the Kuhn-Tucker method depends on the accuracy of the Hessian approximation 
obtained from the RQP method. Using this initial test set we hope to observe how closely 
the Hessian approximation comes to the exact Hessian an d draw some initial conclusions 
on its importance to the accuracy of the results. 

A measure of the closeness of the Hessian approx imation to the true Hessian can be 
defined using the Frobenius norm as 

£ H = II H - H ap prox ll F (5.1) 

This measure has been used in the past to compare the convergence of different variable 
metric updates (Dennis and Schnable 1983). 

For test problem 1 the true Hessian of the Lagran ;ian is 

H _r 2 64 o i 

H “L 0 2.6 J 

From the RQOPT program we obtained the following He ssian approximation with the BFS 
update 

„ _r 1.50017 -.540310 I 
h BFS-L-.540310 2.34388 J 

which gave us a £hbfs = 1-396. 

Using the SRI update from the same starting point, we obtained the following 
Hessian approximation 

„ _r 2.63194 -.002930 1 
H SRl-[. 002930 2.61374 J 

with gives a £hsri = 0.0164. This represents a large improvement in the closeness of the 
Hessian approximation. However, even though the Hest ian approximation for the SRI 
update is much better than the Hessian approximation for the BFS update the problems 
were solved in the same number of iterations (and functk n/constraint evaluations) of 
RQOPT. 

The results given above were obtained with a val re of 8 =1.1. The 8 parameter 
controls the size of the active set during the course of an i >ptimization; a large 8 will cause 
more near active constraints to be considered as part of the active set, a small value of 8 will 
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allow only truly active constraints to be considered ^ of ^ Hessian 

resolved obtaining the following Hessian approximanons 

r2.329 .3227 "| 

Hbfs=|_.3227 2.264 J 

r2.6394 .00093 1 

Hsri-L. 00093 2 5990 -1 ^ 4 were obtai ie d. These improved 

,he values of chbfs - 1 0 60 “ „ was able ,o ideality the correct active set of 

Hessian approxtan^t «n RQ0PT requhed the same number of 

Another implementation issue that the 

Hessian approximation obtained from the opt ”““ d ^ effect of allowing 

reoptimizations performed to esdmate ^ V »•« » <«““ 3 

Hessian updates during the reopumt hied The Hessian approximation that was 

problem 1 was estimated with this option ena . )xima uon that was obtained 

used at the star, of the sensitivity angsts *£ » " ™ we obtained the 
with die BPS update and 5= U. After estimating the sens,, v ty, 

following Hessian approximation 

r2.63975 .00013 1 

H B FS=L .00013 2.59957 J improving the Hessian 

with £hbfs = 0.00053. This indicates that there ,s a posstbuty for unpro 

approximation if we allow abating during the sensiuvity analyse 

Tests for problem 2 were also performed, whose me Hessian of the Ugrangian ,s 
given by 


r 5 1 1 

n= i 5 1 
Ll 1 5 




Ll 1 5 J . . . 0 htain the following value 

Using die starting point provided in j he BFS opdate 

of the Hessian apptoximauon (from RQOPT) when w 

r 4 4976 1.9976 0.49763- 
„ J 1 9976 2.9976 1.9976 

BFS -Lo. 49763 1.9976 4 . 497 * J e ot [he following value of the 

. , „ _ a 000 When we use the SRI update 

with ehbfs ~ 3 AJUU ‘ 

Hessian approximation 
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r4.5 2.0 0.5 1 
H S ri= 2.0 3.0 2.0 
L0.5 2.0 4.5 J 

£ Hsri = 3 0. If we allow the Hessian matrix to be updated while estimating the sensitivity 
of pi with a Ap = 0.0001 we obtain 


Hbfs 


r4.9448 1.0070 1.1749 n 
: 1.0070 4.9964 0.9665 
.1.1749 0.9665 4.3990 J 


Hsri= 


r5 i i 
1 5 1 
Ll 1 5 


with ehbfs = 0.6541 and with e HS Ri = 0.00. This represents a significant improvement 
of the Hessian approximations, particularly when the SRI update is used. 


If we calculate the sensitivity of p2 and use the SRI update we also obtain exact 
convergence of the Hessian approximation. However if we use the BFS update we do not 
obtain exact convergence but an improvement similar to that of the first problem is 
achieved. 


For Test problem 3 the Hessian of the Lagrangian is 


H= 


[-12 0 0 0 I 
0 8 0 0 
0 0 10 0 
. 0004 . 


If we use the starting point that was provided in the problem description, the approximation 
to the Hessian of the Lagrangian (form RQOPT) using the BFS update is 


r 9.785 -0.4657 -2.502 -0.9879 

„ -0.4657 7.7556 -0.5500 0.0174 

HbfS= -2.502 -0.5500 4.062 0.7464 

.-0.9879 0.0174 0.7464 2.1579 

with ehbfs = 7-76. we use the SR1 u P date t0 solve the P roblem 111611 we obtain the 
following approximation to the Hessian 

11.9744 -0.02493 -0.04194 0.02212 
-0 02493 7.9792 -0.03037 0.01095 
Hsri= -0.04194 -0.03037 9.9679 0.00222 
.0.02212 0.01095 0.00222 4.01399 

with a ehsri = 0.130. This represents a major improvement in the Frobenius norm. 


For Test problem 4 the Hessian of the Lagrangian is 


(“6.72 -4.0 
-4.0 9.4006 


H=J 


- 2.0 


6.4 

L-2.0 


- 1.2 

- 6.2 

6.4 


-2.0 6.4 

- 1.2 - 6.2 
4.4 -1.2 

-1.2 9.3418 
-2.0 -4.0 


-2.0 I 
6.4 
- 2.0 
-4.0 
6.26884 
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using the starting point that was defined in the appendix, KQOPT with the BFS update 
yields the following Hessian approximation 



r6.280 

-3.963 

-1.417 

6.458 

-1.924' 


-3.963 

8.052 

-0.561 

-6.247 

5.775 

Hbfs= 

-1.417 

-0.561 

1.548 

-0.932 

-1.449 

6.458 

-6.247 

-0.932 

9.3465 

-4.024 


L- 1.924 

5.775 

-1.449 

-4.024 

5.974 - 


with Ehbfs = 3.6226. 


When we attempted to use the SRI update, the He isian approximation became 
nearly singular after 5 iterations and the Hessian approximation was automatically reset to 
the identity matrix by RQOPT. RQOPT delivered the following Hessian approximation 


["4.661 1 -4.8848 0.1290 5.2854 -3.5129 
-4.8848 7.5178 -.1721 -7.0517 4.7140 
H S ri= 0.1290 -0.1721 1.0046 .1862 -.1245 
5.2854 -7.0517 0.1862 8.6328 -5.0996 
L- 3.5329 4.7140 -.1245 -5.0996 4.4094 


with £hsri = 9.307. The inaccuracy of this Hessian appro ximation results because a total 
of only 7 iterations were needed to solve the problem, and a reset occurred after the fifth 
iteration. Therefore, only 2 iterations could be used to build the Hessian approximation. 
In the near future we will investigate why the Hessian approximation became nearly 
singular after the 5 th iteration. 


A summary of the results of this section are presented in the Table 5.1 where eo 
represents the error between the true Hessian and the identity matrix used at the outset of 
the optimization. Using the BFS update we see that we w ere not able to converge to the 
exact Hessian but the inaccuracies do not seem to be large . As mentioned before, this may 
be due to RQOPT not using exact line searches which the BFS method assumes. Allowing 
updating of the Hessian approximation during the sensitiv ity analysis seems to improve the 
estimate of the Hessian of the Lagrangian. 


Using the SRI update we were able to obtain better estimates of the Hessian of the 
Lagrangian for both problem 1 and 3. For problem 2 the Hessian of the Lagrangian that 
was produced by the SRI update had converged in a projected or reduced sense. The 
inaccuracies in problem 4 are due to a near singular point which is discussed above. 


Problem 

£0 

EBFS 

esRl 

EBFS with updating 

1 

2.291 

1.396 

0.061 

0.0005 

2 

7.348 

3.0000 

3.0 

0.654 

3 

15.93 

7.76 

0.130 

5.186 

4 

23.276 

3.6226 

9.307 

1.698 


Table 5.1 A comparison of the Frobenius norms c f the Hessian approximations 
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5.3. RESULTS 


This section presents a comparison of the sensitivity derivatives calculated by the 
new method with the known sensitivities for the problems in the initial test set We also 
present a means that can be used to compare the accuracy c-f the sensitivity derivatives 
graphically. 

The measure for accuracy that will be used was alto used by Sandgren (1977). 
Sandgren compared the closeness of the optimum design point generated to known 
optimum point, and the closeness of the value of the knowi optimum objective function 
value to the generated value of the optimum plus a penalty for any violated constraints. 
Sandgren defined the following measures 

lABS[f(x)] for f(x* 


) * 0 
I = 0 


(5.2) 


where f(x*) is the true value of the optimum and f(x) is the value returned by the 
algorithm. The total error is calculated as 

nineq neq 

e t = e f + S<gj>+ X(hi) (5-3) 

j=i i=l 

where < a > = (0, if a > 0 I -a if a < 0). The e t measure is ised because it does not bias any 
constraints. 

The relative error in the x vector is defined as 



in equation 5.4 if xj* is equal to zero then the relative erro in xj is defined as the value of 

XJ. 


We will define the relative error in the gradient (df */dp) of the objective function as 

follows 
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f ° r V * ° 


for df 


(5.5) 


We will define the relative error in 3x*/3p and 3u»/3p in the san .e manner as 0 X and denote 
these values as and es 0 ./3p respectively. Eight digits of iccuracy were matnuuned 

in calculating the relative errors. 


The optimal sensitivities for the test problems were calculated using the Kuhn- 
Tucker method with exact derivatives. Once the optimal sensitivities for the problems were 
known, experiments were conducted using RQSEN on the initial test set Both the forward 
difference and central difference variants of the RQP sensitivity algorithm were tested wit 
large and small values of perturbation for the parameters. For nil cases, RQOPT 1 was 
allowed to perform two iterations to optimize the perturbed problem. However, there were 
some instances where RQOPT required only one iteration to meet the convergence criteria. 
A spreadsheet was used to automate the calculation of the relati ve errors in the derivatives 
using the formulas given above. Summary tables showing the relative errors in the 
calculation of the derivatives will be presented for each of the problems. 


Plots of the optimal sensitivities for large variations in the parameters were also 
created for all of the parameter sensitivities that were studied, rhe interesting plots will be 
included in the appendix of this report. These plots can be used to help asses the 
nonlinearity in the sensitivity derivatives, and to help to unden tand the effect of changes in 

the active set. 

The rest of this section presents tables and figures showing the relative accuracy of 
the sensitivity derivatives. A brief discussion of the results for each problem is offered. 


S ^ 1 Problem 1 

Problem 1 possesses three parameters for study. Sensitivities of the objective 
function, design variable, and Lagrange multipliers for each parameter were estimated 
using the four variations of the basic algorithm. The results are compared against the exact 
sensitivities in Tables 5.2 - 5.7. In most cases, the estimated sensitivities agree with the 
known sensitivities with few exceptions. As might be expected, the central difference 
approximations in all cases provides better estimates than the forward difference 
approximations. No strong conclusions with respect to the choice of Ap can be drawn 
from this problem. For parameter 1 , both sires of Ap provide exact sensitivities. For 


ljhe gradients of the objective function and constraints were calculated i sing 
approximations. 


central and forward difference 



parameter 2, the larger value of Ap provides better results wl ile for parameter 3, the 
smaller value of Ap provides the best results. A review of the sensitivity plot for this 
parameter (figure A.2) shows that the sensitivities for this parameter are nonlinear. 

In conclusion, for this problem, using a central difference approximation with either 
step size for Ap resulted in no significant errors in the sensit vity estimates. 


Kuhn Tucker 
Method 



AP =2% relative error AP-0. 1 % 


relative error 


df/dp 1.0000000 

dxi/dp 0.0000000 
dx2/dp 0.0000000 

e x 


1.00000000 0.00E+00 

0.000000 0.000E+00 

0.000000 0.000E+00 

0.00E+00 


1.00000000 0.00E+00 

0.000000 O.OOOE+OO 

0.000000 0.000E+00 

0.00E+00 


dui/dp -0.2000000 
du2/dp 0.4000000 

£u 


-0.200000 0.000E+00 

0.400000 0.000E+00 

0.00E+00 


-0.200000 0.000E+00 

0.400000 0.000E+00 

0.00E+00 


Table 5.2 Central Difference Approximations to the Parameter Sensitivities for problem 1 
parameter 1 



relative error 

-3.00E-03 

0.000E+00 

0.000E+00 

0.00E+00 

-3.159E-03 

7.670E-03 

8.29E-03 


Table 5.3 Forward Difference Approximations to the Pan meter Sensitivities for problem 1 
parameter 1 
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Kuhn Tucker Central Difference Approximations 


Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp -0.3000000 

-0.30000000 

0.00E+00 

-0.30000000 

0.0OE+00 

dxi/dp 0.1000000 
dx2/dp -0.2000000 

0.10000011 

-0.20000030 

-1.100E-06 

3.000E-06 

0.10000359 

-0.20000320 

-3.590E-05 

-1.600E-05 

e x 


3.20E-06 


3.93E-05 

dui/dp -0.1304000 
d U2 /dp 0.0008000 

-0.13040090 

0.00079973 

-6.902E-06 

3.363E-04 

-0.13040262 

0.00080051 

-2.009E-05 

-6.325E-04 

Eu 


3.36E-04 


6.33E-04 


Table 5.4 Central Difference Approximations to the Parameter Sensitivities for problem 1 
parameter 2 


Kuhn Tucker Forward Difference Approximations 


Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp -0.3000000 

-0.30000000 

0.00E+00 

-0.30000000 

0.00E+00 

dxi/dp 0.1000000 
dx2/dp -0.2000000 

0.10004811 

-0.20017630 

-4.811E-04 

-8.815E-04 

0.10000681 

-0.20000665 

-6.810E-05 

-3.325E-05 

Ex 


1.00E-03 


7.58E-05 

duj/dp -0.1304000 
duVdp 0.0008000 

-0.13061831 

0.00114555 

-1.674E-03 

-4.319E-01 

-0.12851547 

0.00999551 

1.445E-02 

-1.149E+01 

Eu 


4.32E-01 


1.15E+01 


Table 5.5 Forward Difference Approximations to the Parameter Sensitivities for problem 1 
parameter 2 
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Kuhn Tucker 
Method 

Central Difference Approximations 
AP =2% relative error t P=0. 1 % 

relative error 

df/dp 

-0.4000000 

-0.40000000 

0.00E+00 

-0.40000000 

O.OOE+OO 

dxi/dp 

dx^/dp 

e x 

1.2000000 

0.6000000 

1.19995460 

0.60007142 

3.783E-05 

-5.952E-05 

7.05E-05 

1.19999990 

0.60000018 

1.000E-06 

-3.000E-07 

1.04E-06 

dui/dp 

du2/dp 

£n 

0.4048000 

-0.5896000 

0.40501856 

-0.58972201 

-5.399E-04 

-2.069E-04 

5.78E-04 

(.40480055 

-(.58960030 

-1.359E-06 

-5.088E-07 

1.45E-06 


Table 5.6 Central Difference Approximations to the Parameter Sensitivities for problem 1 
parameter 3 


Kuhn Tucker 
Method 

Forward Difference Approximations 


AP =2% 

relative error 

AP=0.1% 

relative error 

df/dp -0.4000000 

-0.36400000 

9.00E-02 

0.39820000 

9.00E-02 

dxi/dp 1.2000000 
dx2/dp 0.6000000 

e x 

1.20029290 

0.59484082 

-2.441E-04 

8.599E-03 

8.60E-03 

1.20001670 

0.59973857 

-1.392E-05 

4.357E-04 

4.36E-04 

dui/dp 0.4048000 
duVdp -0.5896000 

eu 

0.39514199 

-0.57797893 

2.386E-02 

1.971E-02 

3.09E-02 

0.40367242 

-0.59207284 

2.786E-03 

-4.194E-03 

5.03E-03 

Table 5.7 Forward Difference Approximations to the Parame er Sensitivities for problem 1 


parameter 3 


5.3.2 Problem 2 

Tables 5.8 - 5. 1 1 present results obtained for the two parameters of problem 2. 
There is some significant disagreement between the known ai d estimated sensitivities for 
parameter 1 in each of the variations tested. The inaccuracies seem to occur due to the 
Hessian approximation not converging. The cause of this is 1 tkely due to RQOPT not 
using exact line searches, which the BFS variable metric update assumes. In this case, 
possibly the SRI update would produce better results. The n '.suits for parameter 2 are 
excellent for all variations. 

The major conclusion to draw from this problem is th it the convergence of the 
Hessian approximation can be a critical factor in the success if the new method. 
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Kuhn Tucker Central Difference Add oximations 


Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp -7.0000000 

-7.00000000 

-1.43E-1 1 

-7.00000000 

-1.43E-11 

dxi/dp 2.0888889 

2.1593080 

-3.371E-02 

2.1598250 

-3.396E-02 

d X2 /dp -2.1666667 

-2.1459360 

9.568E-03 

-2.1459879 

9.544E-03 

dx3/dp -0.9166667 
e x 

-1.0143604 

-1.066E-01 
1.12E-0 1 

-1.0139152 

-1.061E-01 

1.12E-01 

dui/dp -4.6666667 
£u 

-4.5219200 

3.102E-02 

3.10E-02 

-4.5311680 

2.904E-02 

2.90E-02 

dg2/dp -6.00000000 
e e 

-6.1756450 

-2.927E-02 

2.93E-02 

-6.1738965 

-2.898E-02 
2.90E-02 8 


Table 5.8 Central Difference Approximations to the Parameter Sensitivities for problem 2 
parameter 1 


Kuhn Tucker 

Forward Difference Approximations 


Method 

AP=2.% 

relative error 

AP=0. 1 % 

relative error 

df/dp -7.0000000 

-7.00000000 

-1.43E-11 

-7.00000000 

-1.43E-11 

dxi/dp 2.0888889 
dx2/dp -2.1666667 
dxydp -0.9166667 

2.2087420 

-2.1161414 

-1.1367750 

-5.738E-02 

2.332E-02 

-2.401E-01 

2.2379818 

-2.1100841 

-1.1301356 

-5.738E-02 

2.332E-02 

-2.401E-01 

£x 


2.4 8 E -01 


2.48E-01 

dui/dp -4.6666667 

-4.4230730 

5.220E-02 

-4.3803787 

5.220E-02 

£u 


5.22E-02 


5.22E-02 

dg2/dp -6.00000000 

-6.43386710 

- 7.231E-02 

-6.3725934 

-7.231E-02 

e g 


7.23E-02 


7.23E-02 


Table 5.9 Forward Difference Approximations to the Par. meter Sensitivities for problem 2 
parameter 1 
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Kuhn Tucker 
AP =2% 

C pnrral Difference ADDroximations 

Method 

relative error 

AP=0. 1 % 

relative error 

df/dp 12.00000000 

12.00000000 

0.00E+00 

12.00000000 

0.00E+00 

dxj/dp 0.33333333 
dx2/dp 0.33333333 
dx3/dp 0.33333333 

Ex 

0.3333333 

0.3333333 

0.3333333 

1.000E-08 

1.000E-08 

1.000E-08 

1.73E-08 

0.3333333 

0.3333333 

0.3333333 

1.000E-08 

1.000E-08 

1.000E-08 

1.73E-08 

dui/dp 2.33333333 
Eu 

2.3333333 

1.429E-08 

1.43E-08 

2.3333335 

-7.143E-08 

7.14E-08 

dg2/dp 2.00000000 
e g 

2.0000000 

0.000E+00 

0.00E+00 

2.0000000 

O.OOOE+OO 

0.00E+00 

Table 5.10 Central Difference Approximations to the Parameter Sensitivities for problem 2 
parameter 2 

Kuhn Tucker 
Method 

Forward Difference A nnroximations 


AP=2.% 

relative error 

AP=0. 1 % 

relative error 

df/dp 12.00000000 

12.00000000 

0.00E+00 

12.00000000 

0.00E+00 

dxi/dp 0.33333333 
dx^/dp 0.33333333 
dx-j/dp 0.33333333 

e x 

0.3333333 

0.3333333 

0.3333333 

1.000E-08 

1.000E-08 

1.000E-08 

1.73E-08 

0.3333333 

0.3333336 

0.3333333 

1.000E-08 

1.000E-08 

1.000E-08 

1.73E-08 

dui/dp 2.33333333 
Eu 

2.3333333 

1.429E-08 

1.43E-08 

2.3333335 

1.429E-08 

1.43E-08 

dg2/dp 2.00000000 
e g 

2.0000000 

0.000E+00 

0.00E+00 

2.0000000 

0.000E+00 

0.00E+00 


Table 5.1 1 Forward Difference Approximations to the Parameter Sensitivities for problem 2 
parameter 2 

5 3.3 Problem 3 

Problem 3 possesses 3 parameters for study, and he results are presented in Tables 
5.12 - 5.17. Again, the central difference approximation . produce the better results, with 
the small step perturbation better than the large step perturbation for the first two 
parameters. The accuracy of the estimates are good in comparison with the known 

sensitivities. 
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Kuhn Tucker 

Central Difference Approximations 


Method 

AP =2% 

relative error 

AP=0.1% 

relative error 

df/dp -8.0000000 

-8.0000001 

-8.75E-09 

-8.0000001 

-8.75E-09 

dxi/dp -1.1471984 
dx2/dp -0.3342830 
dxydp -0.2279022 
dxVdp -3.5403608 

-1.1551824 

-0.3360517 

-0.2189132 

-3.5312020 

-6.960E-03 

-5.291E-03 

3.944E-02 

2.587E-03 

-1.1485459 
-0.3355475 
-0.2265911 
- 1.5390270 

-1.175E-03 

-3.783E-03 

5.753E-03 

3.767E-04 

£x 


4.05E-02 


6.99E-03 

dui/dp -67.3428300 
du 3 /dp 55.4605800 

-67.5423370 

55.6363930 

-2.963E-03 

-3.170E-03 

-6 7 . 3399900 
58.5155000 

4.217E-05 

-9.903E-04 

£u 


4.34E-03 


9.91E-04 

dg2/dp -1.6600260 

-1.6579776 

1.234E-03 

-1.6595000 

3.169E-04 

e g 


1.23E-03 


3.17E-04 


Table 5.12 Central Difference Approximations to the Parameter Sensitivities for problem 3 
parameter 1 


Kuhn Tucker 

Forward Difference Approximations 


Method 

AP =2% 

relative error 

AP=0.1% 

relative error 

df/dp -8.0000000 

-8.0000001 

-8.75E-09 

-3.0000001 

-8.75E-09 

dxi/dp -1.1471984 
d X2 /dp -0.3342830 
d X3 /dp -0.2279022 
dxVdp -3.5403608 

-1.0434709 

-0.3298696 

-0.2847166 

-3.5084682 

9.042E-02 

1.320E-02 

-2.493E-01 

9.008E-03 

-1.1427923 

-3.3378248 

-0.2292447 

-3.5376210 

3.841E-03 

-1.060E-02 

-5.891E-03 

7.739E-04 

e x 


2.66E-01 


1.27E-02 

dui/dp -67.3428300 
du 3 /dp 55.4605800 

-62.2695100 

52.1090250 

7.534E-02 

6.043E-02 

-66.8998810 

55.1460460 

6.578E-03 

5.671E-03 

£u 


9.66E-02 


8.68E-03 

dg2/dp -1.6600260 

-1.6647592 

-2.851E-03 

-1.6589239 

6.639E-04 

E g 


2.85E-03 


6.64E-04 


Table 5.13 Forward Difference Approximations to the Parameter Sensitivities for problem 3 
parameter 1 
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Kuhn Tucker Central Difference Apj roximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.0000000 

0.0000000 

0.00E+00 

0.0000000 

0.00E+0O 

dxi/dp 

dx2/dp 

dxydp 

dxVdp 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

0.000E+00 

0.000E+00 

0.000E+00 

0.000E+00 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

O.OOOE+OO 

0.000E+00 

O.OOOE+OO 

O.OOOE+OO 

ex 



0.00E+00 


0.00E+00 

dui/dp 

du3/dp 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

O.OOOE+OO 

0.000E+00 

0.0000000 

0.0000000 

O.OOOE+OO 

O.OOOE+OO 

e u 



0.00E+00 


0.00E+00 

dg2/dp 

1.00000000 

1.0000000 

0.000E+00 

1.0000000 

O.OOOE+OO 

e g 



0.00E+00 


0.00E+00 


Table 5.14 Central Difference Approximations to the Parameter Sensitivities for problem 3 
parameter 2 


Kuhn Tucker Forward Difference Approximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.0000000 

-2.200E-13 

2.200E- 13 

-1.660E-12 

1.660E-12 

dxi/dp 

dx^/dp 

dx^dp 

dx^dp 

0.0000000 

0.0000000 

0.0000000 

0.0000000 

4.150E-05 

-3.300E-05 

-2.020E-05 

-3.090E-05 

-4. 1 50E-05 
3.300E-05 
2.020E-05 
3.090E-05 

8.316E-04 
-6.601E-04 
-4.050E-04 
-6.180E -04 

-8.316E-04 

6.601E-04 

4.050E-04 

6.180E-04 

Ex 



6.46E-05 


1.29E-03 

dui/dp 

du3/dp 

0.0000000 

0.0000000 

8.840E-03 

-1.270E-02 

-8.840E-03 

1.270E-02 

1.769E-01 

-2.553E-01 

-1.769E-01 

2.553E-01 

e u 



1.5 5 E -02 


3.11E-01 

dg2/dp 

1.00000000 

1.0000100 

-1.000E-05 

1.0020100 

2.010E-03 

e g 



1.00E-05 


2.01 E -03 


Table 5.15 Forward Difference Approximations to the Parameter Sensitivities for problem 3 
parameter 2 
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Kuhn Tucker Central Difference App roximations 


Method 

AP =2% 

relative eiror 

AP =0.1% 

relative error 

df/dp -10.0000000 

-9.99999986 

1.4E-08 

-9.99999986 

1.4E-08 

dxi/dp 1.3057920 
dx2/dp 0.5460588 
dx3/dp -1.0541310 
dxVdp -2.3741690 

Ex 

1.3055639 

0.5465549 

-1.0523005 

-2.3704118 

1.747E-04 

-9.085E-04 

1.737E-03 

1.583E-03 

2J3E-03 

1.3079005 

0.5481210 

-1.0520338 

-2.3720509 

-1.615E-03 

-3.777E-03 

1.990E-03 

8.921E-04 

4.65E-03 

dui/dp 55.4605800 
du 3 /dp -56.5052230 

Eu 

55.4057600 

-56.4808140 

9.884E-04 

4.320E-04 

1.08E-03 

55.3658650 

-56.9128710 

1.708E-03 

-7.214E-03 

7.41E-03 

dg2/dp 0.67758780 

0.6766357 

1.405E-03 

1.4 IE -03 

0.f.767562 

1.227E-03 

1.23E-03 


Table 5.16 Central Difference Approximations to the Parameter Sensitivities for problem 3 


parameter 3 


Kuhn Tucker Forward Difference Appr QM lanQM 


Method 

> 

II 

to 

cR 

relative error 

AP=0.1% 

relative error 

df/dp -10.0000000 

-9.99999986 

1.4E-08 

-9.99999986 

1.4E-08 

dxi/dp 1.3057920 
dx^dp 0.5460588 
dx3/dp -1.0541310 
dx 4 dp -2.3741690 

£x 

1.3726972 

0.5453465 

-0.9886692 

-2.3442141 

-5.124E-02 

1.304E-03 

6.210E-02 

1.262E-02 

8.15E-02 

1.3151920 
0 5489030 
-1 0457860 
-2 3672333 

-7.199E-03 

-5.209E-03 

7.916E-03 

2.921E-03 

1.23E-02 

dui/dp 55.4605800 
du 3 /dp -56.5052230 

Eu 

56.8153210 

-56.9128710 

-2.443E-02 

-7.214E-03 

2.55E-02 

55 5101970 
-56 7627370 

-8.946E-04 

-4.557E-03 

4.64E-03 

d g2 /dp 0.67758780 
e g 

0.6668760 

1.581E-02 

1.58E-02 

0.6757969 

2.643E-03 

2.64E-03 


Table 5.17 Forward Difference Approximations to the Parameter Sensitivities for problem 3 
parameter 3 


5 .3.4 Problem 4 

The final problem in the test set contains 4 parameters, the estimated sensitivities are 
presented in Tables 5.18 - 5.23. As before, the central difference approximations produce 
the best results, and for this problem, the small step perturbation performs best. Excellent 
agreement was achieved for all the parameters in this problem 
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Kuhn Tucker Central Difference Approximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.517404073 

0.51740407 

0.00E+00 

0.51740407 

0.00E+00 

dxi/dp 

-0.40000000 

-0.4000000 

0.000E+00 

-0.400000 

0.000E+00 

dx2/dp 

0.09729967 

0.0972994 

2.960E-06 

0.097300 

-3.710E-06 

dxydp 

-0.20000000 

-0.2000000 

O.OOOE+OO 

-0.200000 

0.000E+00 

cbcj/dp 

0.28602513 

0.2860244 

2.447E-06 

0.286026 

-3.077E-06 

dxs/dp 

Ex 

-0.06773997 

-0.0677393 

9.994E-06 

1.07E-05 

-0.067741 

-1.247E-05 

1.34E-05 

du3/dp 

0.40709971 

0.4071024 

-6.681E-06 

0.407102 

-6.731E-06 

dus/dp 

-0.05662417 

-0.0566248 

-1.038E-05 

-0.056625 

-1.038E-05 

du6/dp 

0.34730309 

0.3473052 

-6.191E-06 

0.347305 

-6.219E-06 

du9/dp 

6u 

0.01908492 

0.0190855 

-3.280E-05 

3.56E-05 

0.019086 

-3.233E-05 

3.52E-05 


Table 5.18 Central Difference Approximations to the Pan meter Sensitivities for problem 4 
parameter 1 


Kuhn Tucker Forward Difference Approximations 


Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 0.517404073 

0.51740407 

0.00E+00 

0.51740407 

0.00E+00 

dxi/dp -0.40000000 
dxz/dp 0.09729967 
dx 3 /dp -0.20000000 
dxVdp 0.28602513 
dx 5 /dp -0.06773997 

-0.4000000 

0.0973210 

-0.2000000 

0.2860770 

-0.0677900 

0.000E+00 

-2.193E-04 

0.000E+00 

-1.814E-04 

-7.388E-04 

-0.400000 

0.097289 

-0.200000 

0.286000 

-0.067716 

0.000E+00 

1.053E-04 

0.000E+00 

8.713E-05 

3.548E-04 

e x 


7.92E-04 


3.80E-04 

du 3 /dp 0.40709971 
dus/dp -0.05662417 
d U6 /dp 0.34730309 
du 9 /dp 0.01908492 

0.4073857 

-0.0565732 

0.3475276 

0.0190319 

-7.025E-04 

8.997E-04 

-6.464E-04 

2.780E-03 

0.407116 

-0.056623 

0.347315 

0.019083 

-4.073E-05 

2.283E-05 

-3.481E-05 

9.583E-05 



3.07E-03 


1.12E-04 


Table 5.19 Forward Difference Approximations to the Pa rameter Sensitivities for problem 4 
parameter 1 
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Kuhn Tucker Central Difference Approximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.306110869 

0.30611087 

0.00E +00 

0.3061 1087 

0WE+00 

dxi/dp 

0.00000000 

0.0000000 

0.000E+00 

0.000000 

0.000E+0O 

dx2/dp 

-0.14711868 

-0.1471186 

6.797E-07 

-0.147119 

4.758E-07 

dx3/dp 

0.00000000 

0.0000000 

0.000E+00 

0.000000 

O.OOOE+OO 

dxVdp 

-0.04916518 

-0.0491649 

4.841E-06 

-0.049165 

3.620E-06 

dxs/dp 

0.09817962 

0.0981794 

2.332E-06 

0.098179 

1.742E-06 




5.42E-06 


4.05E-06 

du3/dp 

-0.05662417 

-0.0566248 

-1.063E-05 

-0.056625 

-1.038E-05 

dus/dp 

0.05051545 

0.0505156 

-3.662E-06 

0.050516 

-5.424E-06 

du6/dp 

-0.06156394 

-0.0615644 

-6.903E-06 

-0.061564 

-7.115E-06 

du9/dp 

0.00240162 

0.0024015 

3.935E-05 

0.002402 

-2.373E-06 

Eu 



4.15E-05 


1.39E-05 


Table 5.20 Central Difference Approximations to the Para neter Sensitivities for problem 4 
parameter 2 


Kuhn Tucker Forward Difference Approximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.306110869 

0.30611087 

0WE+00 

0.30611087 

O.OOE+OO 

dxi/dp 

0.00000000 

0.0000000 

0.000E+00 

0.000000 

O.OOOE+OO 

dx^dp 

-0.14711868 

-0.1470723 

3.150E-04 

-0.147117 

1.196E-05 

dx^dp 

0.00000000 

0.0000000 

0.000E+00 

0.000000 

O.OOOE+OO 

dxVdp 

-0.04916518 

-0.0490525 

2.292E-03 

-0.049161 

8.675E-05 

dxs/dp 

0.09817962 

0.0980709 

1.107E-03 

0.098176 

4.190E-05 

Ex 



2.5 6 E -03 


9.71E-05 

du3/dp 

-0.05662417 

-0.0570971 

-8.352E-03 

-0.056648 

-4.279E-04 

dus/dp 

0.05051545 

0.0510986 

-1.154E-02 

0.050545 

-5.817E-04 

du6/dp 

-0.06156394 

-0.0620483 

-7.868E-03 

-0.061589 

-4.015E-04 

du9/dp 

0.00240162 

0.0024686 

-2.790E-02 

0.002405 

-1.404E-03 

Eu 



3.23E-02 


I.63E-03 


Table 5.21 Forward Difference Approximations to the Par.imeter Sensitivities for problem 4 
parameter 2 
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Kuhn Tucker 
Method 


df/dp 1.183954566 


dxi/dp 

dx2/dp 

dx3/dp 

dxVdp 

dxs/dp 


- 0.20000000 

0.09408231 

-0.35000000 

0.22881747 

0.02931310 


duj/dp 

dus/dp 

dus/dp 

dus/dp 


0.34730308 

-0.06156394 

0.72069515 

0.15964688 


rvntral nifferenceApprowmati ons. 

AP =2% 

relative error 

AP- 0.1% 

1.18395457 

0.00E+00 

1.1 1395457 

-0.2000000 

0.0940821 

-0.3500000 

0.2288169 

0.0293137 

O.OOOE+OO 

2.689E-06 

0.000E+00 

2.666E-06 

-2.030E-05 

-0.200000 

0.C94082 

-0.350000 

0.228817 

0.029314 


2.06E-05 


0.3473054 

-0.0615644 

0.7206968 

0.1596474 

-6.709E-06 

-8.089E-06 

-2.248E-06 

-3.382E-06 

0. 147305 
-0.061564 
0.720697 
0.159647 


1.13E-05 


relative error 

0.00E+00 

0.000E+00 

2.317E-06 

0.000E+00 

2.273E-06 

-1.723E-05 

1.75E-05 

-6.277E-06 

-7.115E-06 

-2.456E-06 

-3.570E-06 

1.04E-05 


Table 5.22 Central Difference Approximations io the Paramete, 
parameter 3 


Sensitivities for problem 4 


Kuhn Tucker 
Method 


df/dp 1.183954566 


dxi/dp 

dx2/dp 

dxs/dp 

dxrj/dp 

dxs/dp 


- 0.20000000 

0.09408231 

-0.35000000 

0.22881747 

0.02931310 


e x 


du3/dp 

dus/dp 

dus/dp 

dus/dp 


0.34730308 

-0.06156394 

0.72069515 

0.15964688 


fni — if" 1 ^fWnce AppmxmMQS 
AP=2% relative error \P=0.1% 


1.183954566 

0.00E+00 

1.18395457 

-0.2000000 

0.0941409 

-0.3500000 

0.2289599 

0.0291757 

0.000E+00 

-6.226E-04 

O.OOOE+OO 

-6.226E-04 

4.687E-03 

0.200000 

0.094083 

0.350000 

0.228818 

0.029313 


4.77E-03 


0.3485794 

-0.0614102 

0.7229895 

0.1595192 

-3.675E-03 

2.497E-03 

-3.184E-03 

8.000E-04 

0.347369 

-0.061557 

0.720811 

0.159641 


5.52E-03 



relative error 

O.OOE+OO 

0.000E+00 

-2.073E-06 

0.000E+00 

-2.098E-06 

1.559E-05 

1.59E-05 


-1.895E-04 

1.153E-04 

-1.611E-04 

3.620E-05 

2.77E-04 




Table 5.23 Forward Difference Approximations 
parameter 3 


to the Parameter Sensitivities for problem 4 
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Kuhn T ucker Central Difference Ape roximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.010389619 

0.01038962 

0.00E+00 

0.01038962 

O.OOE+OO 

dxi/dp 

0.00000000 

0.0000000 

0.000E+00 

0.000000 

O.OOOE+OO 

dx2/dp 

-0.03975934 

-0.0397594 

-8.803E-07 

-0.039759 

-1.082E-06 

dx3/dp 

0.00000000 

0.0000000 

0.000E+00 

0.000000 

0.000E+00 

dxVdp 

0.07614087 

0.0761408 

1.1 16E-06 

0.076141 

1.366E-06 

dxs/dp 

Ex 

0.15499104 

0.1549911 

-5.162E-07 

1.51E-06 

0.154991 

-6.452E-07 

1.86E-06 

du3/dp 

0.01908492 

0.0190855 

-2.950E-05 

0.019086 

-3.238E-05 

dus/dp 

0.00240162 

0.0024016 

-7.745E-06 

0.002402 

-2.498E-06 

du6/dp 

0.15964687 

0.1596474 

-3.382E-06 

0.159647 

-3.633E-06 

du9/dp 

eu 

0.08386033 

0.0838607 

-4.031E-06 

3.09E-05 

0.083861 

-4.042E-06 

3.29E-05 


Table 5.24 Central Difference Approximations to the Par; meter Sensitivities for problem 4 
parameter 4 



Kuhn Tucker 

Forward Difference Approximations 



Method 

AP =2% 

relative error 

AP=0. 1 % 

relative error 

df/dp 

0.010389619 

0.01038962 

0.00E+00 

0.01038962 

0.00E+00 

dxi/dp 

dx2/dp 

dxydp 

dxVdp 

dxs/dp 

0.00000000 

-0.03975934 

0.00000000 

0.07614087 

0.15499104 

0.0000000 

-0.0397539 

0.0000000 

0.0761542 

0.1549782 

0.000E+00 

1.374E-04 

0.000E+00 

-1.744E-04 

8.265E-05 

0.000000 

-0.039759 

0.000000 

0.076143 

0.154989 

O.OOOE+OO 

1.856E-05 

0.000E+00 

-2.359E-05 

1.123E-05 

e x 



2.37E-04 


3.20E-05 

du3/dp 

dus/dp 

du6/dp 

du^dp 

0.01908492 

0.00240162 

0.15964687 

0.08386033 

0.0193233 

0.0024360 

0.1599642 

0.0841484 

-1.249E-02 

-1.432E-02 

-1.988E-03 

-3.435E-03 

0.019097 

0.002403 

0.159663 

0.083875 

-6.561E-04 

-7.328E-04 

-1.032E-04 

-1.755E-04 

£u 



1.94E-02 


1 .00 E -03 


Table 5.25 Forward Difference Approximations to the Parameter Sensitivities for problem 4 parameter 4 


5.3.5 A Graphical Comparison of the Accuracy Delivered bv Several Different 

Methods 

This section presents a graphical method that can be used to assess the accuracy of 
the parameter sensitivity derivatives calculated by variou: methods. When we plot a bar 
chart of the calculated sensitivity derivatives versus the true sensitivity derivatives, we can 
identify which components of the gradient are least accur ate. A brief discussion of the 
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significance of each graph will be provided. 


Figure 5.1 and 5.2 present a graphical comparison of th< dx/dpi and du/dpi (as 
calculated by various methods*) for problem 3. In addition to results discussed in section 
5.3.3, we calculated sensitivities by the modified Kuhn-Tucker method. The Hessian 
approximations are provided in section 5.2. 

We can see in figure 5. 1 that for large values of Ap for f orward differencing we did 
not obtain results as accurate as those delivered by the other var iants of the RQP sensitivity 
algorithm. If we examine the bars for the KT w/ BFS Hess we see that there is some 
discrepancy in the calculated values of the sensitivity derivative s. We also can see that 
using the Hessian approximation delivered by using the SRI update we were able to obtam 
results better than those obtained by using the BFS approximate on. 



Value of Sensitivity Derivative 


KT 

RQP cd AP=2% 
RQPcd AP=0.1% 
RQP fd AP=2% 
RQPfd AP=0.1% 
KT w/ BFS Hess 
KT w/SRl Hess 


Figure 5.1 A comparison of the accuracy of dx/dp for 3>roblem 3 parameter 1 


We can see in figure 5.2 that all but the KT w/ BFS Hess method produces good 
estimates of the sensitivity derivatives. The discrepancy in du/dp is due to the Hessian 

2 Kp 2 pmblen, wiih cenual difference 

rqp cu »— *•> 

RQP fd AP=2% iSj 2 total' rel« gSbeil I tnbletn wiih forward difference 

RQP fd nP= 0 ..% 1 - , RQP algortthm^^^ 2 it^^tons'tosobrererttuhot problem widi forward difference 
armroximations Ap=0.1% of the nominal value 

KT w/ BFS Hess - The Kuhn Tucker sensitivity equations were solved using the approximate Hessian 
Slivered hv ROOPT when using the BFS update 
KT w/ SRI Hess - The Kuhn Tucker sensitivity equations were solved using the approximate Hessi 
delivered by RQOPT when using the SRI update 
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approximation not converging. 
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Value of Sensitivity derivative 


Figure 5.2 A comparison of the accuracy of du/dp for problem 3 parameter 1 


Figures 5.3 and 5.4 present comparisons of the accuracy of the sensitivity 
derivatives (calculated by various methods 2 ). We can sec that all methods are in good 
agreement with the exception of the calculation of the K1 with the approximate Hessian. 
Again the Hessian approximation did not converge fully for this problem. 
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Value of Sensitivity Derivatives 
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RQP cd AP =2% 
RQPcd AP =0.1% 
RQPfd AP=2% 
RQPfd AP =0.1% 

KT w/ appx BFS Hess 


Figure 5.3 A comparison of dx/dp for Problem 4 as 


calculated by various methods 
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Value of Sensitivity Derivatives 
Hgure 5.4 A comparison of various methods for calculation of 3u/dp for Problem 4 

5 A rONH X TSTONS FPOM TF.ST RESULTS. 

In summary, we feel that this initial testing has shown hat the RQP based method 
can produce reliable estimates of parameter about 

this initial testing. 

We saw in seclion 5, .ha, the Hessian app— 

approximation to be updated dunng the sensiuvity an ys - ng as we calculated 

u, • ,h. t «t set we observed the Hessian approximation improving as 
problems m the test set we obs Ressian arpr0 ximation did not converge 

the parameter sensitivities. This imp . oro ; ecte d sense) then a good 

during the solution of the original problem ” ^ sensitivity analysis. Once we 

estimate of the Hessian approximation can ui witch from a more 

nave a good approximation of .he forwald differonce 

expensive cential differencing approximation <o a less «P«“ condusion is also 

approximation to die sensitivity dertvanve approxma . Kuhn-Tucker sensitivity 

encouraging because the RQP sensitive algorithm approacnes the Kuhn Tuc 

algorithmas tite approximation improve (as was demonstrated in section 3.2). 

.fdreHessianappro— 

T HCSSi "irHeLCCever, Its need to be developed to check if tite Hessian 

updating it while using tite RQP sensitivity algorithm can , ause tire Hesstan » P P 
toconverge. When this happens tite Kuhn-Tucker sensithmy ^ 

the Hessian approximation. Again, the issue to be resolved before dusopuon 


investigated is how to test a Hessian approximation for convergence. 

Our experiments with the SRI update were very encouraging. Section 5.2 
demonstrated how the SRI update was able to deliver a more accurate estimate of the 
Hessian approximation than the BFS update. We also observed exact convergence of the 
Hessian approximation with the SRI update if we allowed th 3 Hessian approximation to be 
updated during the sensitivity analysis. Near convergence ol the Hessian approximation 
for problem 1 and 3 was also obtained when the SRI update was used. 

Even though our experiments with the SRI update were encouraging the SRI 
update has some serious drawbacks that will have to be investigated further before the 
method can be recommended. Unlike the BFS update the SRI update is not self-correcting 
(Ip 1987). As we saw in problem 4, the SRI update delivered an Hessian approximation 
that was singular, and the RQP method performs best if the Hessian approximation is 
positive definite. 

In section 5.3, we observed that the method was able to approximate sensitivity 
derivatives. We saw that using central differencing yielded more accurate estimates of the 
sensitivity derivatives than using forward differencing appn >ximations. We did not 
observe any major sensitivity to Ap when we were using the central differencing option. 

We also observed that the forward differencing approximation has more trouble 
evaluating sensitivity derivatives, and the perturbation step size Ap had more effect when 
the functions being approximated are nonlinear in the parameter. 

Finally, we observed that when the sensitivity derivative is small, the RQP 
sensitivity algorithm can sometimes have a difficult time finding the exact value. This was 
demonstrated in the calculation of 3u2/3p for problem 1 anr. also in the calculation of the 
sensitivity derivatives for p2 in problem 3. 
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6 . Changes in the Active Set 


This chapter describes our findings regarding changes in the active set of 
constraints. The first section describes the four cases of changes in the active set of 
constraints. The second section provides sample plots illu strating the predicted behavior 
from section 1. Next, the third section discusses several proposed solutions for dealing 
with changes in the active set. The final section provides tome procedures that can be used 
to generate problems where there are changes in the active set. 


6. 1 . CASES OF ACTIVE SET CHANGE 

When the active set of constraints change some of he sensitivity derivatives may be 
discontinuous. It is also possible that at the point where the active set changes, the 
gradients of the constraints may become linearly dependent. When this happens the 
optimization problem may become very difficult to solve. The four cases for changes in the 
active set are; 

1. A constraint enters the active set (constraint gra iients are linearly independent). 

2. A constraint leaves the active set. 

3. A constraint enters the active set (constraint gradients are linearly dependent). 

4. A constraint enters the active set and causes the re to be no feasible solution 

Case 1 and Case 2 are complementary 1 . This can (test be demonstrated by an 
example. Consider an optimization problem whose solution changes as a parameter p, is 
varied and the active set changes as p increases. Assume the original problem is optimized 
for a p=p° and and a sensitivity analysis is performed. The value of p is then increased to 
p=p! (where p 1 > pO) and the problem is reoptimized. At p=pl, constraint g 3 will enter the 
active set if p is increased any further, but the gradients of the constraints remain linearly 
independent as p is increased. Now p is increased to p=p- (where p2 > pi > pO) and at this 
point constraint g 3 has been added to the active set A Case 1 active set change has 
occurred. To see the complementarity, assume we begin with p=p2 and conduct a 
sensitivity analysis for decreasing values of p. Now when p = pi, for any value of p less 
than pi constraint g 3 will leave the active set. Thus, when p reaches pO, constraint g 3 will 
have left the active set and we have had a Case 2 change ir the active set for p going from 
p 2 to pO. 

The behavior of the sensitivity derivatives for Case s 1 and 2 is characterized as 

1 This means that active set change algorithm that are developed for ase 1 changes in the active set may 
also be used for case 2 changes in the active set 
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follows. There is a discontinuity in the rate of change of the objective function with respect 
to p (d 2 f*/dp 2 is discontinuous), there can be a discontinu ty in the rate of change of the 
some of the design variables (9x*/5p is discontinuous), at d the rate of change of the 
Lagrange multipliers can be discontinuous (i.e. 5u*/5p for a constraint leaving the active 
set will go from some nonzero value to zero, this may cau ;e a change in the rate of change 
of the other Lagrange multipliers as well). The Hessian of the Lagrangian is continuous 
with respect to variations in p. For these problems there will be directional derivatives of 
dx*/9p and du*/dp that can be used to make estimates of how the design variables and 
Lagrange multipliers will change. Second order extrapolations of the behavior of the 
objective function can also be made using d 2 f*/dp2 in a di ectional sense. 

The characteristic of Case 3 is a discontinuity in th<: Lagrange multiplier estimate 
(u* is discontinuous). The discontinuity in the Lagrange multiplier causes a discontinuity 
in df*/dp and also causes the Hessian of the Lagrangian t< » be discontinuous. Since the 
active set changes, there will also be a discontinuity in (d> */dp). At the point where the 
constraints become linearly dependent it will not be possible to solve the Kuhn-Tucker 
sensitivity equations because they will become singular. For Case 3, there may be an 
exchange of constraints in the active set (i.e. the new cons raint may replace one of the 
constraints that is already in the active set as p moves through the point). 

The main characteristic of Case 4 is that there only exists a directional derivative 
away from the point where the path terminates. In order t< < calculate the directional 
derivative for Case 4 we need to be able to find the proper active set to follow when we 
leave the degenerate point. Case 4 changes in the active sot will be common when the user 
overconstrains the design, i.e. sets the performance specif ications to high to be physically 
meet by the given design configuration. 


6.2. DEMONSTRATION OF CASES 

This section describes the behavior that we observ id in the initial test set for 
problems that had changes in the active set. This section dso discusses two test problems 
that are not in the initial test set that are used to demonstra e the behavior of Case 3 changes 
in the active set. 

The problems in the initial test set only had Case 1 and Case 2 changes in the active 
set. Indications of where the active set changes were sh< >wn on the plots of the optimum 
sensitivity provided in the appendix. 

If we examine the plots that are presented in the a rpendix (figures A.1-A.8, A. 10, 
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and A.l 1) we can see the discontinuity in 9x*/9p and 9u*/9p when the active set changes. 
In problems 3 and 4, because there is only a small discontinuity in 9x/9p, it is difficult to 
see in the graphs provided (figures A.5,6,7,10,1 1) in the appendix. However for 
problems 3 and 4 it is easy to see the discontinuity in 9u/<>p, the slope of the Lagrange 
multipliers, when the active set of constraints changes. 

If we examine figures A. 4 and A. 1 1, we see that r ot all of the components of 
dx*/dp are discontinuous when the active set changes. For problem 2 in figure A.4 we do 
not see a discontinuity in 9x2*/9p, and for problem 4 in figure A. 10 we do not see a 
discontinuity in 9xj*/9p or 9x3*/9p, when the active set ciianges. This behavior is partially 
due to the fact that we were studying right hand side perti rbations for these problems. The 
symmetry in problem 2 can be used to explain its behavior. In problem 4, the gradient of 
constraint gg had zero's in the first and third locations, thi s we expect that perturbations of 
constraint gg should have no effect on 9xi*/9p or 9x3*/9p. 

The discontinuity in d 2 f*/dp 2 is difficult to see, when the active set changes for 
most of the problems in the initial test set. However for test problem 2, parameter 1, 

(figure A.3) we can see the discontinuity in d 2 f*/dp 2 . In the plot of the optimum value of 
the objective function versus pi we can see that for value ; of pi < 1.055, d 2 f*/dp 2 is less 
than zero (a region of negative curvature) and for values of pi > 1.055, d 2 f*/dp 2 is greater 
than zero (a region of positive curvature). 

In problem 2 (p 2 ), problem 3 (p1.p2.p3) and problem 4 (p4) we studied the effect of 
perturbing the right hand sides of the inequality constraints. In these problems we studied 
a range of perturbations from where the constraint was inactive to where the constraint is in 
the active set. When the constraint was inactive perturbations in the parameters had no 
effect on the optimum. When the constraints entered the active set for these problems they 
caused d 2 f*/dp 2 to go from zero to some nonzero number The discontinuity of d 2 f*/dp 2 
can be seen in figures A.4, A.5, A.6, A.7, and A.l 1. 

The following example will demonstrate how the optimum behaves for Case 3, 
when the gradients of the constraints become linearly dependent upon a change in the active 
set of constraints. 

minimize: f = xi2 + (P - 1)2 (6. i) 

subject to: gi = 3 xi + 2 P - 10 > 0 (6.2) 

g2 = 2 xi + 3 P - 10 > 0 (6.3) 

when P = 2, the minimum f* = 5 occurs at xi* = 2 with both constraints active and the 

gradients of the constraints linearly dependent. The Lagra lge multipliers and in the family 
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ui,u2 e {3 ui + 2 U 2 = 4, ui > 0, U2 > 0} ' 

For this problem df*/dp, ax*/9p and 3u/a P do not exist. T1 ese derivatives do exist in 
both the positive and negative directions and are indicated b / dx/dp+ for increasing values 
of p and 9x/dp- for decreasing values of p. 


Figure 6.1 presents the sensitivity information for p roblem 6. 1-6.3. Figure 6.1 (a) 
and (b) represents the first order predictions of the new values of the Lagrange multipliers 
for this problem. For this problem the linear predictions agree with the optimum Lagrange 
multipliers. There is a discontinuity at Ap = 0.0, therefore mere are only be directional 
derivatives for these values. Figure 6.1 (c) represents linear predictions of the new value 
of the objective function. Notice again that there is a discontinuity in the slope of the 
prediction. Thus df*/dp does not exist for this value of p (however df*/dp + and df*/dp- 
will exist in a directional sense where the superscript + or - ndicate the direction of change 
in p). Figure 6.1 (d) represents the predicted location of the optimum. Notice the 
discontinuity in the slope at Ap = 0.0, the true location of the optimum agrees with the 

linear prediction for this problem. 



Figure 6.1 Plot of the optimum sensitivity predict ions for problem 6.27-29 

We have also examined a test problem used by Fi icco and Ghaemi (1982). This 
test problem has a linear dependence in the constraints gradients. The problem is to design 
a corrugated bulkhead for an oil tanker. The objective fur ction is to minimize the weight of 
the bulkhead with constraints placed on allowable bending stress, moment of inertia, 
corrosion, and minimum gage thickness of the material. The problem has six design 
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variables, 18 design parameters, 17 constraints, a non i J 

linear consuls. We did no. apply any scaling to dns problem, even .hough 

constraints for this problem are poorly scaled. 

A sensitivity analysts was performed for variations in th. parser, LT, front 

0*1* » 1* <*«* - RQ0 “^ "tion value, design 

convergence of the objective function value, pied J ^ analysis was 

variables, and die Lagrange multipliers were also u . ns were performed 

performed before the RQS y variables, Lagrange 

by hand. Therefore we did not record all of the values oi 

multipliers, and constraints. 

Figure 6.2(a) is a plot of the sensitivity of the optimal objective function for the 
bulkhead problem versus the parameter LT which was we have 

482.8 we notice a discontinuity in the slope of e ^ is a Case 2 change in 

encountered a Case 3 change in die active set. A. LT - 476 . 2 . 

the active set. 

,, oresents a pi 0 , of the value of constrain. 13 (g,3> versus the value of LT. 
A, LT- 482.8, the value of this constnun, goes to zero and the contain, enters the active 

set. 

Figure 6 3 presents a plot of die value of the Ugrang, multiplier for constiain, 12 

s LT Notice that there is a discontinuity a. LT - 482.8 this is due to a linear 
versus LT. Notice tnai uic a«TT moves from one side of 

s::it,:rwt rat:: ^ se, * .^ts w„ h gI3 a 

constraint on the lower bound of x 6 . 

Figure 6.2(d) presents a plot of die optimal value of design ^ 

482.8, the lower bound constraint is active. However, w e plots of the 

changes. 
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g(13) u(12) and u(13) 



Figure 6.3 Plot of The Lagrange Multipliers foi the Bulkhead problem 



Figure 6.4 Plot of gi 3 versus Lt for the Bulkhead problem 


For this problem there is a second change in the active set at LT = 476.25 when g 12 
leaves the active set. As LT moves through 476.25 there s not the same type of 
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discontinuity as we saw for Case 3. This represents Casn 1, a constraint leaving the active 
set. We can see some of the characteristics of this case i i the plots. In figure 6.2, the 
slope of the objective function is continuous at LT = Alt. 25. In figure 6.3, the value of the 
Lagrange multiplier of constraint gj 2 goes to zero signifying the constraint has left the 
active set. Finally, in figure 6.2 we can see that 5xj/3p, 9x3/9p, and dx^/dp are 
discontinuous at LT = 476.25. 

A Case 4 change in the active set can be illustrated by the following simple one 
variable problem. 

minimize f(x) = x 

Subject to: gi = x - pi > 0 
g2 = x - 2 > 0 

When pi = 1 the optimum solution to this problem is f*=T.O, x*=1.0, and constraint gi is 
active. When pi is increased to pi=2.0 constraint g 2 enters the active set. However, if pi 
is increased beyond 2.0 there is no feasible solution for this problem. 

An engineering example of a case 4 change in the active set can be illustrated by the 
following example. Find the optimum air plane to fly a g iven mission, where the design 
variables may be the size of the wings, engine size, cruise altitude, etc. The design 
parameters might be; total cargo weight, runway length, air temperature at take-off, etc. 
Assume a constraint on take-off distance is applied. 

Take-off distance < Available Runway Length 
and a parameter sensitivity study of total cargo weight is performed. As the total cargo 
weight is increased, the take-off distance of the airplane increases, and the design variables 
of the airplane may also change. If the total cargo weight is increased beyond a certain limit 
it may become impossible for the aircraft to take-off at the given runway length. A 
parameter sensitivity analysis for values of the total cargo weight for values greater than 
this limit are meaningless because there will be no solution to the problem since the plane 
cannot take-off. Thus, for this problem when the runway length constraint enters the active 
set there will be a case 4 change in the active set causing no feasible solution to exist. 


6.3. PROPOSED SOI .1 JTTONS 

This section will discuss some algorithms that car be used to obtain more accurate 
estimates of the parameter sensitivity after the active set of constraints has changed. The 
first subsection describes algorithms that can be used to calculate sensitivity derivatives 
when there are Case 1 or Case 2 change in the active set. The second subsection presents 
an example demonstrating how the algorithms presented n the previous section work. The 
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final subsection discusses how sensitivity derivatives can be calc elated for Case 3 changes 
in the active set. 

derivatives. 

Directional derivatives can be approximated as follows 


dfW) .Limit ( ■ f (p° + Ap) • f*(P°r (6.5) 

dp Ap-^0 + ^ Ap 

Here dfW)/dp indicates the rate of change of a function when a parameter is perturlwl 
Here df* ( P V P im fa # analysis is degenerate, we 

in the positive direction. „ lh - ROP sensiti'ity algorithm in the same 

can approximate directional denvanves usmg the RQP senstn tty 

way as we approximate other derivatives. 


6 3 1 Doling wit h 1 an( * — 

We begin this section by presenting the following example problem 

minimize: f = (M + ' “' 2 

subject to: gi =xi -P-° 

g2 = 6 - 2 xi - x 2 > 0 


( 6 . 6 ) 

(6.7) 

( 6 . 8 ) 


. r £ c f * _ 4 v* = (1,2) A sensitivity analysis 
When d = pO = 1 as shown in figure 6.5, f - 4, x t eventually 

P .. „ j = (10) Notice, that as we increase p, eventually 

indicates that as p is mcrea , »P ^ ( active ^ Increasing p 

constraint g2 will become acuve intersection of the two 

from (x*,pl) will result in the optimum pot*" and dten 
constraints. The deflection algorithm ts based on ftndmg a cons™" 
calculating new values of 3**/3p and 3u*/3p along the new acuve . 
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Figure 6.5 Deflection of the Sensitivity Se arch Direction 


The proposed algorithm for dealing 


with constraints entering or leaving the active 


set is given by the following steps. 


1. 

2 . 


ermine dx*/3p at the base point 
culate Apl for intersection of the constraint 
a mnstraint to enter/leave the active set 


by finding the minimum Ap that 


. i min ( gj ) ; g Active set 

* ■ 


Ap 


1 = 


mm 


( u; 




j e Active set 


3. Calculate dx* Vdp with the active set updated to 

step 2. 

4. Calculate Ax and Au by 


(6.9) 

(6.10) 


reflect the change indicated by 


Ax= *37 Apl + 

r)ll* 

Au ^Apl + 


3x*i . 0 

3ll*l A n O 

5r Ap2 


( 6 . 11 ) 

( 6 . 12 ) 


where 
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Apl = min ( p - pO, 0) 
Ap2 = max ( p - pi, 0) 


(6.13) 

(6.14) 


5. Calculate the new value of the objective function by 

fnew = f(x*,p0) + Ap |^+ j^-Ax (6.15) 

The first step in the above algorithm is to perform a sensitivity analysis to determine 
which direction the optimum will move in. In the second step, we attempt to identify the 
value of p for which the active set will change using equa tions 6.9 and 6.10. A new search 
direction is calculated in step three which includes the changed active set In step four, a 
formula is provided which defines the proper value for Ax for before the constraint is 
encountered and after the constraint becomes active (inactive). This result is used in step 5 
to predict a new value for the optimum at the new point. 

In step three of the deflection algorithm, it is necessary to find the deflected search 
direction at the point where the new constraint enters the t ctive set. This can be done by 
adding the newly violated constraint (or removing the constraint that is leaving the active 
set) to the Kuhn-Tucker sensitivity equations and then solving for dxVdp and 3uVdp. 

This computation can be performed efficiently by using matrix updating techniques as 
described in (Diewart 1984). In order to obtain a more accurate estimate of the 
sensitivities, the gradients of the active constraints can be re-evaluated at x = xl, where x* 
is the predicted intersection of the new constraint. 


When using the deflection algorithm, a check shot Id be made to assure that the 
constraint gradients are not linearly dependent. If the constraint gradients are linearly 
dependent then this procedure cannot be used to solve for :he optimum sensitivities because 
there is no solution to the Kuhn-Tucker sensitivity equations. The procedure to be used 
when the constraint gradients are linearly dependent will !* discussed in section 6.3.3. 


The estimate obtained in step 5 is a linear extrapolation. A better estimate of the 
new value of the optimum can be made using a quadratic e xtrapolation. The formula for a 
quadratic extrapolation is (McKeown 1980, Fiacco 1983, Barthelemy and Sobieski 1983) 


,* £/ jj,-. . df* 1 . d 2 f* . 

f new - f(x ) + A Pi gpr + y A Pi A Pi 


(6.16) 


where 


d 2 f* 

a 2 L 

- . i 

a 2 L 1 

T 3x* 

5uT 

3vT dn 

dpi 2 

" a P i 2 ! 

_dpidx*_ 

^pT 

'cJpT 5pi 



(6.17) 
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Using (6.16) to predict the new value of the objective function as p varies yields a 
more accurate estimate of the new optimal objective function until the active set changes. 
We should note that equation 6.17 is in terms of 9g/dp and dn/dp, and can be rewritten in 
the following form in terms of 8x*/8p 


d 2 f* 8 2 L 
dpi 2 9pi 2 


2 


■ a 2 L ]T dx*_ ^ 8x*T 8 2 L ax* 

.a Pi ax*. api + api ax * 2 ^pi 


(6.18) 


If we substitute equation 6.18 into equation 6.16 and write the equation in terms of 
Ap and Ax we will obtain 


i r a 2 L a 2 L T 

f*new = fl st order + 2 dpi 2 + dpidx 


a 2 L 




] * 


19 ) 


where the first order approximation comes from equation 6 15. This form has the 
advantage of being able to predict the new value of the objective function at a new value of 
x and p. When the deflection algorithm is used to calculate the new location of the 
optimum then equation 6.19 can be used to make a second order estimate of the new value 

of the objective function. 


As mentioned in section 6.1, a constraint leaving the active set (Case 2) is 
complementary to a constraint entering the active set (Case l). It may be possible to predict 
df*/dp when a constraint leaves the active set without calcu lating dx^/dp and duVap along 
the new active set This could be beneficial for instances that only require estimates of 
f*(p). To predict the behavior of the optimum when a constraint leaves the active set we 
can use directional derivatives. Using the following formula (Fiacco 1983) 


df* af ^ ^h. 
dP = 3p fl '^P 


nineq ^ 

£, > 


( 6 . 20 ) 


and the fact when an inequality constraint leaves the active set its Lagrange multiplier goes 
to zero, we can obtain an estimate of df*/dp at the point where the change in the active set 
takes place, this is done by obtaining an estimate of 8u/dp from a sensitivity analysis at the 
base point using it be used to estimate when a constraint w 11 leave the active set. A 
prediction of the df*/dp when the active set changes can then be made from the following 

formula 


df* 3f idhi 

dp" = 3p i?i ‘3P 



( 6 . 21 ) 


where 


v i =v 


| *+^ L Ap 1 


( 6 . 22 ) 
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(6.23) 


u ! = Ui * + |r Apl 

When these values of ul and vl are used we obtain a new estimate for the value of the 
df*/dp that will be valid when the value of p is increased or decreased past p 1 . 

The deflection algorithm can be used to predict when a second change in the active 
set will take place. A prediction of when a second constraint may be added to or removed 
from the active set can be made by making a linear approximation using the following 

formulas 


gj =gj + 





e Active set 


(6.24) 



e Active set 


(6.25) 


and 


Ap r= 


= min 


9g; 3gj T 3x* 1 


j e Active set 


(6.26) 


( l A 


Apj =min! 


u 


v 


ia 




j € Active set 


(6.27) 


where Ap 2 predicts the value of pi that will cause the second constraint to enter (using eq. 

6.26) or leave (using eq. 6.27) the active set, gj is a prediction of the value of the j* 

constraint, and uj* is the predicted value of the Lagrange multiplier when the active set 
changes. The predicted value for Ap2 is calculated as a linear approximation of the value 
of the constraint in the dxl/dp direction. If the constraints are interrelated (i.e. are 
evaluated as a set) then it may be possible to use a more ac curate estimate of dg/dx in 
formula 6.26, by predicting its new value by using the formula 


(P 1 ) = (p0) + Ap i 


(6.28) 


When the constraints are interrelated 3 2 g/dx*dp can be ev iluated when dV x L/dp is 
evaluated. 
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If the calculation of the gradient of a particular constraint is not expensive then the 
value of the new gradient may be used in formula 6.26. 

The deflection algorithm and the variants that wer ; introduced in this section will be 
effected by any nonlinearity that is present in the problem particularly nonlinearity in the 
constraints. It should be emphasized that these sensitivities are only estimates of how the 
sensitivity will change. 

6.3.2 An Example 

This section presents an example problem to demonstrate how the deflection 
algorithm of section 6.2. 1 performs. Example problem (5.6-8) will be used. Figure 6.5 
shows the solution for p = 1.0, this example assumes thai p is increased to p = 2.5. The 
exact value of the new optimum at p = 2.5 is f* = 13.25, X* = (2.5,1). 

For p = p0 = 1.0, the initial search direction for an increasing p is dx*/9p = (-1,0) 
as shown in figure 6.5. Step 2 of our algorithm determines that constraint g 2 enters the 
active set when p=p 1=2.0. For values of p greater than pi the location of the optimum is 
along the intersection of constraints gi and g 2 - The new search direction along constraints 
gl and g 2 is determined from step three to be, 

"5p” = C" 1 ’ 2 ) (6.29) 

For p = p2 = 2.5, which is greater than pi, the estimated Ax is composed of the 
sum of two vectors: one vector from (x*,p0) to (x*,pl) plus the vector from (x*,pl) to 
(x*,p2). Thus by equation 6.11 for a Ap=1.5 

Ax = (1.5,-1) (6.30) 

Thus the estimate of the new location of the optimum for i>=2.5 becomes 

x *new est = (2.5,1) (6.31) 

Which is the true location of the new optimum for this problem. Without using the 
deflection algorithm we obtain the following estimate (by equation 2.9) of the optimum 

x *new est = (2.5,2) (6.32) 


By equation 2.7 and 2.8, the new value of the objective function will be 

f|CTi st = 10 (6.33) 

Using equation 6.10, which takes into account the change in the active set, we get 
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fKTDl st = 10 (6.34) 

Equation 6.16, which gives a second order estimate without taking into account the change 
in the active set, we obtain 

fKT2 nd = 12.25 (6.35) 

Equation 6.19, which provides a second order estimate but does take into account the 
change in the active set, we obtain 

fKTD2 nd = 13.25 (6.36) 

Using the deflection algorithm the predicted location of the optimum was in exact 
agreement with the true value. The linear predictions of r* new gave the same value because 
there was no component of the gradient of f in the X2 direction. Using equation 6.16 for 
the second order estimate provided a better estimate of the new value of the objective 
function, but when equation 6. 19 was used the exact value of the objective function was 
obtained. We should probably not expect results this got<l in more general optimization 
applications. However, we can expect better predictions of the location of the new 
optimum and the value of the objective function for small changes in the parameter when 
the active set changes. 

6.3.3 Dealing with Case 3 

Case 3 is the most difficult case to deal with for clianges in the active set because 
when the active set changes, the Lagrange multipliers will be discontinuous and predicting 
the new active set as the parameter moves through the de generate point is very difficult. 

To deal with Case 3, we propose avoiding the singular point by reoptimizing the 
problem for values of p that are on either side of the singular point. Performing a 
sensitivity analysis at both points, use these sensitivities in a directional sense for p moving 
away from the singular point. 

Reoptimizing the problem on either side of the sir gular point may be a difficult 
problem if the point is nearly singular. This may cause the algorithm being used in the 
reoptimization to fail or converge very slowly (Powell 1985, Bartholomew-Biggs 1986). 

When using an algorithm such as RQOPT that updates an approximation of the 
Hessian of the Lagrangian, it may be unclear which values of the Lagrange multipliers to 
use if the Lagrange multipliers are not unique at the solution. 
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6 4 OF.NEP ATTON OF TE ST PROBLEMS 


This section will describe a procedure that can be used to generate test problems that 
possess changes in the active set. The first secdon describes 'he generaaon of test 
problems with Case 1 and Case 2 changes in the active set. Ihe second section describes 
the generation of test problems with Case 3 changes in the ac ive set. 

a 4 i Test Problem s for Case 1 and Case 2 Chan ge? in the Active Ssl 

Several test problems that exist in the literature posse >s Case 1 and 2 changes in the 
active set (See Schtnit and Chang 1984, Vanderplaats and Y. .shtda 1985, Vanderplaats and 
Cai 1987) To generate new problems, active set changes for these two cases can be 
intreduced into a test problem by adding constraints that are not active at the optimum but 
are violated for small changes in the parameters. The follow ing is an oudtne of the steps 
that can be used to generate test problems where there are ch inges in tile active set 

Step 1 Generate a NLP Test Problem. One such method would be to use the 
Rosen and Suzuki (1965) procedure. 

Step 2 Vary some problem parameters and find the rath of the optimal solution 
with respect to p, say x»(p) from pO to p 2 . 

Step 3 Choose a value of p' between p° and p 2 as the point where a constraint will 
enter the active set. 

Step 4. Construct a new constraint such that 

gnew(x*(p 1 )iP 1 ) = ° 
gnew( x *(P l - Ap),? 1 - Ap) > 0 

and the gradient of gnewtxVLP 1 ) * lineajl V '"dependent of the gradients 
of the constraints that are already in the acti ve set. 

Step 5. Calculate the path of the new problem fron t p 1 to p 2 

The above algorithm can be used to create test problems for Case 1 when a 
sensitivity analysis is conducted a. p = p° and then p in penurbedtop 2 . Thesametest 
problem can be used as a Case 2 test problem, when the sensitivity analysts is performed at 

p = p 2 and p is moved to p°. 

This procedure is illustrated in figure 6.6 which shows a graph of a two variable 
test problem along with the optimum x*(P»). At the optimum, constraint g, ts active. As 
the p increases from pO to pi the location of the optimum moves from x*(pO) to x (p ), 
and constraint g, remains active. To introduce a change i i the aenve set we can place an 
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additional constraint g2 (shown in figure 6.6 (b)), that in ersects the vector {9x(pO)/dp}Ap, 
where Ap = p2 - pO. The value of p where the path to the optimum intersects the new 
constraint will be denoted as pi. By adding constraint g2, the optimum location x*(p2) 
shown in figure 6.6 (b) will be different than without the constraint. 




Figure 6.6 The creation of Test Problems with changes in the Active Set 


6.4,2 Generation of Test Problems for Case 3 

This section will describe the generation of test problems that have a Case 3 change 
in the active set We first discuss problems that are in the literature. Then we will discuss 
two different algorithms that can be used to generate test problems of this type. 

A survey of the literature revealed several test proolems where the constraint 
gradients are linearly dependent when the active set changes (Case 3). Three such 
problems were found in articles by Bartholomew-Biggs ( 1986), Vanderplaats and Yoshida 
(1985), and Fiacco and Ghaemi (1982). It is suspected that problems discussed by 
Robertson and Gabriele (1987) and Barthelemy and Sobieski (1983) also possess Case 3 
changes in the active set. Powell (1985) has studied the performance of RQP methods 
when the gradients of the constraints are linearly dependent. The test problems that were 
used by Powell can also be modified to be sensitivity test problems. 

We can expect to find other test problems (Case 3 ) when we begin to study more 
engineering test problems. Many engineering optimization problems are fully constrained 
at the optimum. When a new constraint becomes active for a fully constrained problem we 
will either have a Case 3 change in the active set or loose the feasible region(Case 4). 

To generate test problems for Case 3 changes in ti e active set the algorithm 
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described in section 6.4. 1 can be used with the following modification. In step 4 when the 
new constraint is added to the active set it will have to be linearly dependent with the 
constraints in the active set. To accomplish this the gradient of the new constraint at x(p^) 
can be constructed as a linear combination of the gradients o r the other active constraints. It 
may be possible with further work to control which constraint in the active set is replaced 
when the active set changes. 

An alternative, less general procedure that can be used to create test problems with 
linearly dependent constraint gradients is illustrated for a two variable problem in figure 
6.7. Figure 6.7 (a) shows a simple two variable optimization problem. The problem has 
an elliptic objective function and is subject to an equality constraint that changes as the 
parameter p changes. There are also variable bounds presei t. For p = pO the optimum is 
located at the intersection of the equality constraint and the v ariable bound X2max* When p 
is changed to p = pi 111 the optimum is then located at the intersection of the equality 
constraint and both X2max and xi m i n . At this point the gradients of the constraints are 
linearly dependent which causes the linear independence as* umption of the Kuhn-Tucker 
conditions not to hold. When p changes further to p = pi as shown in figure 6.7 (b) the 
optimum is now located at the intersection of the equality constraint and ximin* Thus as p 
moves from pO to pi there will be a change in the active set with the constraint gradients 
being linearly dependent. This procedure can be generalize 1 to more dimensions. The 
discontinuity in df*/dp can be modified by varying the ecce itricity of the ellipses. 



Figure 6.7 The Creation of Test Problems with Linear Dependencies 



7. 


Conclusions and Recommendation: for Future Research 


7.1. CONCLUSIONS 

In this work we have proposed a new method for estimating parameter sensitivity 
based on the Recursive Quadratic Programming method. The new method approximates 
the sensitivities using a differencing formula and can be shown to be equivalent to a 
modified Kuhn-Tucker method. The method appears to lie very competitive with existing 
methods when measured in terms of the number of function evaluations required to 
calculate a parameter sensitivity. It does not require the calculation of second order 
derivatives, but uses the BFS method or SRI method for developing an approximation to 
the Hessian of the Lagrangian. 

The choice of the variable metric update (BFS or SRI) effects the amount of work 
required to solve the perturbed problem, because different updates provide Hessian 
approximations of differing accuracies. Different variable metric updates also effect the 
speed with which the RQP algorithm can solve the problem. 

Initial testing of the algorithm against problems with known sensivities has shown 
that the method can adequately estimate the derivatives. " Tie central difference 
approximation seems to provide the best results, particularly when the Hessian 
approximation is updated during the RQP iterations at the perturbed point. 

Using the RQP method to solve the optimization problem may be beneficial in many 
applications. This is because the RQP method has been s hown to be one of the best 
general purpose algorithms for solving nonlinear programming problems. The use of 
variable metric updates allow the RQP method to solve pioblems where the Sequential 
Linear Programming (SLP) method fails. The RQP method can solve problems with very 
nonlinear constraints, and the RQP method performs the !>est when there are many active 
constraints at the optimum of the problem. 

Chapter 6 has discussed potential problems occur ing when there are changes in the 

set. Chapter 6 also presented some techniques that can be used to deal with these cases. 

We observed discontinuities of the sensitivity derivatives in Chapters 5 and 6 when there 

were changes in the active set of constraints. If we are using sensitivity derivatives to make 

extrapolations, we now know Case 3 will cause the largest discontinuities (i.e. step 

discontinuities in the Lagrange multipliers) in the derivati ves, and Case 4 will cause no 

solution to the proposed constraints to exist. We have also shown that the discontinuities 

in dx*/9p and 3u*/9p, that occur when the active set changes, make our prediction of 

which constraint will enter the active set second very difficult without reoptimizing the 
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problem. 

7.2. RECOMMENDATIONS FOR FUTURE RESEARCH 

The next step in the testing of the new RQP sensitivity algorithm will be to expand 
the test set to include a larger variety of test problems. We will need to include problems 
that have more variables, and also problems that demonstrate the behavior associated with 
Case 3 and Case 4 changes in the active set. We will also teed to expand the test set to 
include engineering test problems for which the Hessian o the Lagrangian is not readily 
available. 

As we observed, the perturbation Ap can effect the accuracy of the derivatives that 
we calculated and we will have to investigate methods to i mprove the choice of Ap. In our 
initial testing we always let RQOPT use two iterations to solve the perturbed problem. 
Further tests are needed investigate the effectiveness of the algorithm when RQOPT is only 
allowed one iteration to solve the perturbed problem. 

We observed that the Hessian of the Lagrangian improved if we allowed the 
approximation to be updated during the reoptimization. More experiments need to be 
conducted to find how to best update the Hessian approximation. We also need to find 
when we should switch from using two iterations to solve the perturbed problem to using 
one iteration to solve the perturbed problem. Further work is also needed to identify when 
the Hessian approximation has converged. If the Hessian approximation converges to the 
true Hessian of the Lagrangian then the Hessian approximation may be used in the Kuhn- 
Tucker sensitivity equations, however 9V x L/dp will still rued to be calculated. Since 
V X L = 0 the calculation of 3V x L/3p may be subject to numerical noise and the RQP 
sensitivity algorithm may perform better than the Kuhn-Ti cker method. 

The initial testing of the SRI update was very encouraging in terms of the 
convergence of the Hessian approximation to the True Hessian. However the SRI update 
proved to be unstable when used for test problem 4. We will need to investigate some 
method that can be used to stabilize the SRI update or find a set of rales that can be used 
that will only allow the SRI update to be used when the uixiated Hessian approximation 
will be stable. 

Using the RQP method to approximate the Hessian of the Lagrangian may be 

improved further if a Hybrid MOM/RQP algorithm is used to solve the original 

optimization problem. A Hybrid MOM/RQP algorithm would use the Method of 

Multipliers (or Augemented Lagrangian) algorithm for the first few stages, (build an 

approximation of the Hessian of the Lagrangian) then swi ch to using the RQP method. 

When the RQP method is used the approximation of the F essian of the Lagrangian from 
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then MOM method can be used as the initial Hessian approximation, this should help the 
RQP method quickly solve the problem, and also obtain a *ood approximation of the 
Hessian of the Lagrangian. 

In summation we have seen that the new RQP sens itivity algorithm can find 
parameter sensitivities. We still need to investigate if this method will be superior to the 
Kuhn - Tucker method for general problems. Section 3.2. 2 demonstrated that there will be 
a trade off with regarding the number of function evaluations required by the Kuhn-Tucker 
method versus the RQP method. We will have to conduct more experiments to study the 
general accuracy that can be expected from the RQP sensit vity algorithm. 

Using first order extrapolations can provide unsati factory results when the 
functions are nonlinear with respect to p. This situation is illustrated in Figure 7.1, which 
shows the effect of variations in p 3 on xi*. On this plot a inear and quadratic extrapolation 
are presented as well as the actual values of the optimum of xi*(p 3 ). It can be seen that the 
linear extrapolation does not provide a good estimate of the new value of xi* when p 3 
changes, however the quadratic approximation provides an accurate estimate of xi*(p 3 ) up 
until the point where the active set changes. 


If good estimates of the second derivatives can be found then more accurate 
estimates of the behavior of the optimum can be made by using second order 
extrapolations. 

There are few available methods to calculate second derivatives of the optimum with 
respect to parameter variation but it is possible to predict <l 2 f*/dpi 2 for some problems. 
However the only published algorithm that was found for calculating d 2 x*/9pj 2 requires 
third order derivatives, which are seldom available for engineering problems. Thus, an 
algorithm based on the central differencing variants of the RQP sensitivity algorithm is 
proposed that can be used to calculate second derivatives. 


When using the central difference approximation e n estimate of the second 
derivatives can be calculated from 


d 2 y(pj) _ y + - 2y + y - 
dpi 2 (Api) 2 


(7.1) 


Where y can represent f*, x*, u*, or g*, and the estimates of f + »‘, x + >\ u + >~, and g + ’' 
calculated in steps 2 and 4 of the RQP sensitivity algorithm are substituted appropriately in 
(7.1). When using a central differencing approximation this option can be effective for 
indicating when curvature is present, but may not be able io accurately predict the true value 
of the second order information. It should be noted that mis procedure may have the most 
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difficulty in predicting 82u(pj)/dpi2 because the RQP algor thm produces more accurate 
estimates of f + >" and x + >' than u + ’". 

Second derivatives might also be useful for predict ng when constraints enter the 
active set. A prediction of when a constraint enters the act ve set can be identified using the 
following model of the behavior of the constraints that are not in the active set 

gnew = g(x*,p) + ^T A Pi + I^ (Api)2 (7 * 2) 



Figure 7.1 A Comparison of Quadratic to Linear Ext apolations for the Sensitivity 

Approimations. 

As a final note to this report we will discuss using the RQP method and the RQP 
sensitivity algorithm in multilevel decomposition algorithms. Because multilevel 
decomposition algorithms solve the same subproblems for different values of system level 
parameters, we expect that as the multilevel decomposition algorithm converges (after the 
proper active set has been identified) that the true Hessian of Lagrangian in the subsystem 
will not change very drastically. If the RQP method is used to solve the subsystem level 
optimizations and perform the sensitivity analysis, we exj*ect that the Hessian 
approximation will converge to the true Hessian of the La grangian as the multilevel 
decomposition converges. If we use the approximation of the Hessian of the Lagrangian 
from the previous subsystem optimization as the initial Hessian approximation for the next 
iteration, a better approximation of the Hessian of the Lagrangian should be obtained when 
we solve the new subsystem problem. 

The above discussion implies that using the RQP nethod in conjunction with the 
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multilevel decomposition method may mean that on the fin t few system level iterations the 
sensitivities calculated at the subsystem level may be inacc irate, but the accuracy of the 
sensitivity derivatives will improve after each iteration at the system level. In our opinion 
this behavior is suitable for use with multilevel decomposition since far from the optimum it 
is often not advantageous (or necessary) to perform exact 1 ne searches or have exact values 
of the gradients. However, as the solution is approached, ve need to obtain more accurate 
derivatives and perform more exact line searches to obtain acceptable convergence. 
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Appendix 1 

Test problems 


This section will present a discussion of the test p roblems that were used in the 
initial testing. Selected plots of the optimum sensitivity ere provided to show how the 
optimum varies when the parameter is perturbed. The plots will also demonstrate how 
changes in the active set effect the optimum sensitivity. 

On each plot that is presented, the base value of ti e parameter will be indicated by a 
vertical line indicating the value at which the sensitivity analysis was performed. Active set 
changes will be indicated by a vertical line and a label indicating the constraint that enters or 
leaves the active set when the parameter is perturbed from the base point. On the plots of 
the optimum objective function vs the parameter, a linear extrapolation will be indicated by 
a line showing the predicted value of the objective function using equation 2.8. 

The plots of the behavior of the optimum can be used as a tool when diagnosing the 
behavior of various algorithms used to predict parameter sensitivities. Nonlinearity in the 
paths f*(x*,p), x*(p), and u*(p) can be seen in these plots, this nonlinearity can be used 
to explain why (or how far) linear extrapolations are valid. The discontinuities in the 
sensitivities that occur when the active set changes can al so be seen and we can use these 
discontinuities to establish regions where the extrapolations are valid. 
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Problem 1: 

Minimize f(x) = (xi - pi) 2 + (x2 - P3) 2 

2 

subject to: gi(x) = 2xi - X£ - P2 ^ 0 

g2(x) = P3 - 0.8xj - 2x2 > 0 

Variable bounds £qJ < x < 


Starting Point for Optimization x° = (0,3) P° = (3,1,3) 
Optimum Point: f(x*(p°)) = 1.25 x*(p°) = ( 2.5,2.0) 

Both constraints are active: u*(p°) = 0.3, 0.4) 


Hessian of the Lagrangian H= 

[2.64 0 I 

L 0 2.6 J 


Sensitivity derivatives 



q 

1—H 

II 

dx* rO.Ol 

Lo.oJ 

q '4 ] varied pi from 1.5 to 4.8 

-r— = -0.3 
dp2 

3x* r-0. 1 1 

L0.2J 

3u* _ f-0. 13041 

3p2 L 0.0008 J 

varied pi from -3.0 to 1.5 

-r—= 0.4 

dp3 

0x* rl .21 
^ ~ L0.6J 

3u* r .4048 1 
’5p3~ L-0.5896 J 

varied pi from 2.2 to 3.95 


Special features: Hessian matrix is available. Quadratic objective function and 

quadratic constraints. Active set changes are introduced for large changes in 
the parameters. Problem is fully constrained at the optimum. Plots of the 
behavior of the optimum are presented in figures A. 1 and A.2. 

Constraint gi leaves the active set when p 1 > 4.5, and constraint g2 leaves 
the active set when pi < 2.0. 

Constraint g2 leaves the active set when J 3 > 3.888, and constraint gi 
leaves the active set when P3 < 2.4322. S nee u*(p3) is nonlinear a linear 
prediction of when the constraint will lea\ e the active set will not be accurate 
for large variations in the parameter. 
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Figure A.l Sensitivity of Problem 1 with respect to P(l) 







Figure A. 2 Sensitivity with respect to p(3) for Problem 1 



p(3) with P(2)=1.0 mu! P(l)=10 





Problem 2: 


Minimize 



[5 1 1 I 


- 5 ‘ 

0.5x1 

1 5 1 

x + x T 

5pi 


Ll 1 5 J 


L 5 J 


subject to: gi(x) = pixi + X 2 + X 3 - p 2 ^ 0 
ffoCxi = xi + 2x9 + 3x9 - 4.7 - 1 


Variable bounds 


r- 10-1 


- 20 l 

-10 

< X < 

20 

L- 10 J 


l20J 


Starting Point for Optimization x° = (1.1, 1.2, 1.3) p° = (1,3) 
Optimum Point: f(x*(p®)) = 25.5 x*(p®) = ( l.C.1.0,1.0) 

gl active at the optimum: u*(p°) = ( 12 . 0 , 0 . 0 ) 


Hessian of the Lagrangian H=| 


5 1 
1 5 
LI 1 


1 I 
1 

5 . 


Sensitivity derivatives 



3x* 

3pT 


2.0833333-1 
-2.1666666 
. -.9166666 . 


[ lk66 ^ 666 ] varied pi from 0.8 to 1.2 


df 

dp2 


12.0 


a * 0.3333331 

S2- = 0.333333 
d P2 L0.3333333. 


r2.333333l 

L 0 J 


varied p 2 from 2.7 to 3.2 


Special features: Hessian matrix is available. Quadratic objective function and 
linear constraints. Parameter 2 is a right ha nd side perturbation. 

A plot of the behavior of the optimum with aspect to pi is presented in 
figure A.3. For pi > 1.04889, constraint g 2 enters the active set. We can 
see a change in sign of d 2 f*/dp^ from a region of positive curvature to a 

region of negative curvature. We can also S'ie the discontinuity of the slope 
of dx*/dp and du*/dp when the active set changes. 

A plot of the behavior of the optimum with respect to p 2 is presented in 
figure A.4. For pi < 0.95, constraint g 2 enters the active set.. Again we 
can see the discontinuity in 5x*/9p2 and 3it*/9p2- 


C ' 
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Design Variables 










p(2) with pOH 






Problem 3 (Common name Rosen and Suzuki Tert Problem): 
Minimize f(x) = x^ + x 2 + 2x^ + x 4 - 5xi - 5x2 - 21 x 3 + 7x4 + 100 

2 2 2 2 

subject to: gi(x) = (-Xj - x 2 - x 3 - x 4 - xj + X2 - <3 + X4 )/8 + pi > 0 
g2(x) = (-Xj - 2 x 2 - X3 - 2x 4 + xj + xO/10 + P2 ^ 0 
g3(x) = (-2xj - x 2 - X3 - 2xi + X2 + X4)/5 + p3 > 0 


Variable bounds 


r-lOi 


[201 

-10 

<x < 

20 

-10 

20 

L-10J 


L20J 


Starting Point for Optimization x° = (0.0,0.0,0.C ,0.0) p0 = (1,1,1) 

Optimum Point: f(x*(p 0 )) = 56.0 x*(p0) = ( 0 . 0 , 1.0, 2.0, -1.0) 

gl and g 3 active: u*(p°) = ( 1 . 0 , 0 . 0 , 10 . 0 ) 

-12 0 0 0 

Hessian of the Lagrangian H= q q q 

.0 0 0 4 


Sensitivity derivatives 

dx* 


df_ 8n 

357 " - 80 


^ =0 ° 35T’ 


r -1.1471984 H 
-0.33428030 
-0.2279022 
L -3.5403608 J 


rO.Oi 

0.0 

0.0 

L0.0J 


3u* 

3pT' 


9u* 
c>P2 ' 


r-67.34 28300n 
0.0 

. 55.4605800 J 


o.o- 

0.0 

.0.0. 


varied pi from 0.86 to 1.12 


varied P 2 from 0.85 to 1.15 


df - inn dx * ■ 

357 - - 10 - 0 357 


1.30579201 


0.5460588 

3u* 

-1.0541310 

3pi ~ _ 

L-2. 3741690.1 



-67.J428300H 

0.0 

55.^605800. 


varied pi from 0.86 to 1.12 


Special features: Hessian matrix is available. Quadratic objective function and 

quadratic constraints. Problem has been used in the other studies of optimal 
design, the parameters that we are studyin g are right hand side 
perturbations. 
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FigureA.5andA.6presentplo.softhevariatioaoftheopdn.umw,* 

respect to p,. The following changes in *e active set are introduced: when 
p, < .855, consnain. g 3 leaves the active set; v hen p, 2 1.057 constratn. 
g2 enters *e active set; when p, > 1.078, cons, rain, g, leaves *e active se . 

The plot shows that P(p,) is nonlinear. It is difficult to see *e 

discontinuity in dx*/3p, bu, we can clearly set *e discontinutty tn 3u Id Pi- 
Figure A.7 presents a close up view of the bef avior of x 3 and X4- We can 
see that x 3 is a nonlinear and a piecewise continuous function o pi. e 
*scontinuities lake place when the active se. changes. Wnh *e resoluuo 
in figure A.5 *is behavior was very difficult to see, however wnh *e 
enlarged view this becomes easy to see. 

A plot of the behavior of *e optimum wi* re tpect to P 2 is present in 

figure A.8. When p 2 S 0.9, constraint g2 enters the active set. In figure 

A 8 we can see a sharp discontinuity in 3u/3f and we also can see *a. large 
errors would be introduced if the value of the objective funcuon was 
extrapolated from *e base point to after *e a stive set changed. 
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Design Variables K**) 


Figure A. 5 Sensitivity with respect to p(l) for Problem 3 
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Lagrange Multipliers 


Figure A. 6 Sensitivity with respect to p(l) for Problem 3 (cont) 



p(l) with p(2)=p(3)=l 



p(l) w th p(2)=p(3)=l 
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Figure A.7 Plots of the Sensitivity of x(3) an i x(4) for problem 3 p(l) 
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Design Variables 


Figure A.8 Sensitivity with respect to p(2) for Problem ~ 



n 1 1 S 


* x(l) 

* x(2) 

♦ x(3) 

♦ x(4) 


0.85 0.90 0.95 1.00 1.05 1.10 

p(2) wilh p(l)=p(3)=! 


<2 0.20 

c 




0.85 0.90 0.95 1.00 1.05 

p(2) wilh p(l)=p(3)=l 




.80 0.85 0.90 0.95 1.00 1.05 

p(2) willi p(l)=P(3)=l 










Problem 4: 


5 5 5 5 , 

Minimize f(x) = Xejxj + I IcijXjXj+ £d;x j 

j=l i=ij=i i=l ! 


subject to: 


5 


gi(x) — XajjXj * bj > 0 
j=l 


i=l,10 


Where the values of ay, bj, Cjj, dj, ej are constants that can be found in (Coville 
1969, Himmelblau 1972, Eason and Fenton 1974, or Sandgren 1977) 


Variable bounds 


r°- 


r 20 i 

0 


20 

0 

< x < 

20 

0 


20 

LoJ 


L20J 


Starting Point for Optimization xO = (0.,0.,0.,0.,1.) p<> = (b3, bs, b6, b9)=(-.25,-4.,-l.,5.) 
Optimum Point: f(x*(pO)) = -3.234866708 

x*(p°) = (0.3, 0.3334676, O.sO.4283 10 1,0.2239649) 
g3> g5» g6i and g9 are active: 

u*(p°) = (0.,0.,0.517404,0.,i).30611 1,1.183954594,0.,0.,0.010390,0.) 


Hessian of the Lagrangian 


T6.72 -4.0 -2.0 6.4 
-4.0 9.4006 -1.2 -6.2 
H= -2.0 -1.2 4 4 -1.2 

6.4 -6.2 -1.2 9.3418 

L-2.0 6.4 -2.0 -4.0 


-2.0 “I 
6.4 
- 2.0 
-4.0 
6.2688J 


Sensitivity derivatives 


df 

dpi 


= .517404 



-0.4 -I 
0.097300 
- 0.2 

0.286025 
L-0. 067740-1 


r 0.407099-1 
3u* _ -0.056624 
^7' 0.347303 
0.0 1 9084. 


varied pi from -0.39 to -.019 
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df 

dp2 


.306111 


8 x * 

<Jp2 


r o.o i 


-0.056624-1 

-0.147119 

8 u * 

.0505155 

- 0.0 

8 p 2 ~ ~ 

-.061564 

-.049166 1 

.0.00241 )16J 

0.09817 - 




varied P2 from -4.5 to -3.5 


df 

dp3 


1.18395 


8 x * 

^P3* 


-0.2 “1 


.34303 I 

0.094082 

8u* 

-0.061564 

-0.35 


0.720695 

0.228817 

. 0.159647 J 

-0.029313- 




varied p3 from -1.2 to -0.8 


iL = .010390 
dp4 


8 x * 


0.0 i 


-0.019085-1 

-0.039759 

8u* 

0.002402 

0.0 

dp4 

0.159647 

0.0761408 

.0.083860J 

L. 15499104- 




varied p 4 from 4.0 to 6.0 


Special features: Hessian mams is available. Cubic objective function and linear 
constraints. Active set changes are introduced for changes m p 4 . The 
parameters being studied are right hand side perturbation. Hus problem has also 

been studied by Fiacco et. al. (1974,1983). 

A plot of the sensitivity of the optimum with respect to pi is presented m igure 
A.9. Plots of the optimum sensitivity with respect to P2 and p 3 are similar to 
those to pi. For variation in the first three parameters, the optimum objective 
function behaves linearly and there are no changes in the active set for the range 
of parameter perturbation that was studied. 

A plot of the sensitivity of the optimum with respect to pa is presented m figures 
A 10 and A.l 1. Constraint g 9 leaves the active set when p4 ^ 4.876. A ter 
constraint g 9 has left the active set further varia lion of p 4 has no effect on the 
optimum. L increastng values of p 4 , ermrs would * introduced if a linear 

extrapolation was used. 
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Lagrange Multipliers Design Variables 


Figure A. 9 Sensitivity with respect tc p(l) for Problem 4 
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Design Variables 



p(4)=b(9) with p(l)=-.25, p(2)=-4, p(3> =-l 
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Appendix 2 


THE ROCRE AND ROSEN SYSTEMS 

This section will discuss the system th it was created for studying parameter 
sensitivities. The first subsection introduces the support program RQCRE that was written 
to simplify the construction of test problems. The second section will discuss the RQSEN 
system. The RQSEN system is a interactive program that acts as a post 
processor/sensitivity analysis module for the RQOPT program. 

RQCRE and RQSEN were set up to have some user friendly features. The 
programs have some user-freindly features however the programs will crash if, the user 
enters numbers in response to a promt for an <dpha data type, and may crash if the user 
enters letters when the system is requesting numbers. 

A2.1 The RQCRE Support System 

The RQCRE program was written to reduce the time required to implement test 
problems. The RQSEN program requires approximately 30 arrays and a complicated main 
program to be written by the user. The RQCF E program automatically dimensions the 
proper arrays and writes the required calling p rograms. Using the RQCRE program 
reduced the time required to implement test programs during our initial testing. 

The RQCRE program is essentially a ■ irogram that writes another program. The 
main features of the RQCRE program are; 

1. The program can be used in an inte ractive mode. 

2. The program writes the main calling program. 

3. The program can write an outline of the function subprogram. 

4. The program can be used to update the problem formulation. 

A structure chart of the RQCRE progiam (in CMS) is presented in Figure A2.1. 

The basic functions of the program modules sre; 

RQCRE.EXEC - this module connects the proper output files to the proper unit 
numbers. 

RQCRE - is a FORTRAN program that can be used interactively to create a problem 
for submission to the RQSEN system. The input to RQCRE can either 
come from a data file or from he user. The output from RQCRE is a data 
file "data" that can be used as in input file for RQSEN, a MAIN FORTRAN 


1)9 



program ready for compilation, \ nd a shell for the function subprogram. 



cn o o o CM CM CM 



(fsub ij 

Figure A2.1 A Structure Chart foi the RQCRE Program. 


=N 

=LFSW 
=NEQ 
=LEQ 
=NINEQ 
=LINEQ 
=NPARM 
15 =MAXIT 

1 =1 START 

1 =IFIN 

50 =LPP 

6 =N0UTFL 

5 =ITR 

1 . e- 4 =PMIN 
l.E-8 =CAPFMN 

0 . 5 =DELTA 

10.0 =R 

.0001 = GAMMA 

1.0E-8 =EPSQP 
5. E-4 =EPSGRD 

2 =IDIF 

0 =NSCALE 

1. E-11 =ZER0M 

1.1 =X(1) 

1.2 =X(2) 

1.3 =X(3) 

-10. =XMIN ( 1 ) 

-10. =XMIN ( 2 ) 

-10. =XMIN (3) 

20.0 =XMAX (1) 

20.0 =XMAX (2) 

20.0 =XMAX (3) 

1.0 =PARM(1) 

3.0 =PARM(2) 

Figure A2.2 A Sample of a Data file (Test Problem 2) for the RQSEN Program. 
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data - the input/output data file that contains the algorithm parameters, starting point 
and initial values of the parameters for RQSEN. 

MAIN - FORTRAN program used as he main calling program when running 
RQSEN. 

FSUBI - a FORTRAN function that th ; user is required to modify by adding the 
definitions of the constraints ar d objective function. 

A sample of the data file for test problem 2 is presented in Figure A2.2. This file 
was written by the RQCRE program. This da? a file is used as an input to the RQSEN 
system to provide the programs with the values of the algorithm parameters, design 
variables, and initial values of the parameters. 

A sample of the a program written by the RQCRE program is provided in Figure 
A2.3. The program represents an implementation for test problem 2 (described in appendix 
1). The only modifications that were made to he program are, the objective function and 
the constraint definitions that were added to the code generated by RQCRE. 

A 2. 2 The RQSEN program 

This section describes the implementation of the RQSEN system. The first topic to 
be discussed is the capabilities of the RQSEN system. Next a description of the 
implementation is provided. The final topic p resented in this section is a sample session 
from the RQSEN system. 

The basic capabilities of the RQSEN s /stem are; 

1 . The program can be used to solve < iptimal design problems. 

2. The program can be used to conduct convergence studies for various versions of 
RQOPT. 

3. The program can be used to calculate parameter sensitivity derivatives. 

4. The program can be used to conduct studies of large variations in the 
parameters. 

5. The program can be used to create aarameter sensitivity plots that can be used, 
for trade off studies. 

The RQSEN system is currently imph mented on the following systems, an IBM 
4341 under the CMS operating system and a microVAX under the VMS operating system. 
The programs are written in FORTRAN 77 ai d implemented in double precision. 

Figure A2.4 presents a structure char for the RQSEN program (CMS 
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PROGRAM RTS 02 


TEST PROBLEM 2 PH. D., CREATED BY TODD J. BELTRACCHI TO EXAMINE 
CONVERGENCE OF THE HESSIAN APPROX AND CHANGES IN THE ACTIVE SET 
THIS HAS A QUADRATIC OBJECTIVE FUNCTION AND LINEAR CONSTRAINTS 

IMPLICIT REAL* 8 (A-H,0-Z) 

DIMENSION X (3) ,XMIN (3) , XMAX (3) , SCALE (3) , H (3, 3) , DELFHG (3, 3) , 

1 FUNCT (3),P(3),U(8) ,V(1) ,XS (3) , Dl'HGS (3, 3) ,FUNCTS (3) , PS (3) , US (8) , 

2 VS (1) , DFDP ( 2 ) , FUNCTP (3, 2) ,DXDP (3,2) ,DUDP (8, 2) ,DVDP (1,2), DFDPE (2) , 

3 DGDP (2,2) 

LOGICAL YESNO, ACTPT ( 8 ) , ACTPTS ( 8 ) 

CHARACTER* 3 YSN 
CHARACTER*? FILENM 
COMMON /OPTDAT/ D(400) 

COMMON /BFSWK/ DD(20) 

COMMON /PMINI/ PMINI (3) 

COMMON /PMAXI/ PMAXI(3) 

COMMON /PARMS/ PARM(2) 

INCLUDE (RQS) 

END 

C ************************************, ******************************** 
FUNCTION FSUBI (X, IEVAL) 

IMPLICIT REAL* 8 (A-H,0-Z) 

DIMENSION X (3) 

COMMON /NFEVAL/ NCE , NFE 
COMMON /PARMS/ P(2) 

IF (IEVAL .GT. 1) GOTO 2 

1 NFE=NFE+1 

C PLACE OBJECTIVE FUNCTION DEFINITION FERE 

FSUBI=2 . 5* (X (1) **2+X (2) **2+X (3) ** 2) +X (1) *X (2) +X (1) *X (3) +X (2) *X(3) 

1 +5. * (X (1) +P (1) *X (2) +X (3) ) 

RETURN 

2 NCE=NCE+1 

C PLACE LINEAR EQUALITY CONSTRAINTS HEFE 
GOTO (10, 20, ) IEVAL 

10 FSUBI=P (1 ) *X (1) +X (2) +X (3) -P (2) 

RETURN 

20 FSUBI=X (1) +2 . *X (2) +3. *X (3) -4 . 7-F (1) 

RETURN 

END 

Figure A2.3 A Sample Program (Test Problem 2) for RQSEN. 
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PARVAR 


1 


WRTPT 

^r 




[ prvrdr 


RQOPTs 

Routines 


Figure 


(fsubij 

A2.4 A Structure Chan for the RQSEN System. 









implementation). A brief explanation of each of the program modules is provided 

RQSEN.EXEC - this module connects the proper files to the proper unit numbers 
and prompts the user for the name of the program to be run. 

MAIN - the main calling program. This module calls RQOPT and RQSEN, this 
module also reads in the starting point and algorithm parameters for RQOPT 
and allows the user to save solution point to a data file for later use, i.e. a 
sensitivity study at a later time . 

FSUBI - the function subprogram tha: defines the objective function and 
constraints. 

data - the data file that contains the algorithm parameters and values of the design 
variables. 

RQOPT - implementation of the RQP method described in Beltracchi and Gabriele 
(1987 b). 

RQSEN - the main driving routine for a sensitivity analysis. 

WRTPT - A utility program for writing the design point, Lagrange multipliers and 
values of the objective function and constraints to a summary file. The 
summary file is in a form which can be read by a plotting program to 
graphically display the sensitivity information. 

solf - a data file used to communicate optimum design points to a plotting program. 

PARTP - a subroutine that uses either forward, central or user supplied routines to 
calculate the partials with respect to the design parameters of the objective 
function and constraints. 

PARSEN - a subroutine that calculates dx*/3p, 3u/3p, df/dp and dg/dp by either 
forward or central differencing. The user can specify the perturbation size 
and the number of iterations that RQOPT uses to solve the perturbed 
problem. 

PARSEN2 - a subroutine similar to P ARSEN but this imnplements modified central 
differencing. 

RQDR - is the routine that controls the execution of RQOPT for reoptimization. 

PARVAR - the subroutine used to co nduct studies of large variations in problem 

parameters, if parameter sensitivity derivatives are available then the location 
of the starting point for the reoptimization is approximated by. 

^initial = X (P) + Ap 

PRVRDR - is the routine that control > the execution of RQOPT for reoptimization. 



The firs. Step is for the user ,o cmate the „ A2 . 3 ). The 

main calling program and die funcnon subprogram '* A2.2) with the 

second step is to set up a data file (similar to die - P— » * can * 

values of die design variables and algorithm parame era. These fust 

performed with the aid of the RQCRE preproc defined, the 

Once die cafiing program, function of thc opi * of the 
user can run die RQSEN system to conduct a study of die sensrnv y 

problem or to study the convergence of die problem. 

can be used to assess the convergen i . the test problem faster than the version 

A 2, shows diattheBFS vision ^ ^^^^muccurateestimateofthe 
!iprimum,^d once^te^gion of die minimum wa‘ located the convergence for the SRI 

update ^ow toe program can rhe 

sensidvity studies is to provide a sample session in Fig u re A2.7. 

IBM version of RQSEN. The next several paragraphs descn 

all user responses are shown in italics. . lvcd The modified SRI 

update is^twedto !^pr^^te tte^wsianof 

first prompt asks the user for a n summary of all optimum points. The 

with rt . 02. This data file can be used to main, « m02 , an 

next prompt is for the The program and data file 

implementation of test problem 2 (see appendix 

were presented in Figure A2.3 and A2.2. of the with t he input 

When die then reads in die dau from file rt sO, 

rn^^ltLlta^msevemlarra^ Thesearmys 
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Value of Design Variable Relative Error 



0 10 20 30 40 50 60 70 80 90 100 

Iteration 


figure A2.5 A Sample Plot comparing the Convergence Characteristics of the BFS and SRI 


U pdates. 



0 10 20 30 40 50 60 


Iteration 


Figure A2.6 A Sample Plot Showing the Convergence of the Design Variables. 
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are saved if the optimum has already been found, a response of "n " for no ,s entered . 

Now the program invokes the RQOPT program to locate the optimum of the problem, 
summary of the output from RQOPT (during the opdmiaadon) is presented^ figure 
A2 7(a-c). After the problem has been solved som: final statistics from RQO “« 
presented along with an approximation to the Hess ,an of the Lagrangian, m both die LDL 
format and the unfavored form. These are shown in Figure A2.7(b & c). 

After the problem has been solved by RQOPT control is returned to the 

preprocessor (see Figure A2.7(d). The preprocessor provides the user with a choice o 
being able to save the optimal point In the example the response was ' y for yes was 
entered, next the user is asked if he wants to save the final point in the same file as die 
initial point, a response of "n " for no was entered. Next the user is asked to supp y a 
new name of the data file to store the point in, a response of rts02s was entered.. = 
data file was then written and the user asked if they wanted the gradients and Hessian 
approximation to be written to the file, a response of "y " was entered The prep~ 
next asks if the user wants to perform a sensitivity analysis, a response of y 

enKred Now control is passed the the RQSEN program (see Figure A2.7(c-e)). The first 
question asked by RQSEN is if the user wants the solution points written to a so ution i e a 
response of "n ’ is entered. The next question asked is for a value of epsilon to calculate 
die partial derivatives of die problem functions. RQSEN then calculates die partial* of the 
objective function and constraints, the derivative of the objective function is then calculated 

by equauon iJO.^ ^ ^ outpgt is fte calculation of the partials of die optimum 

design variables and Lagrange multipliers with respect to the first parameter for the 
problem. Again the user can specify the size of he perturbation of the P^ met “ and * 
number of iteration that RQOPT is allowed to use for solving the perturbed problem, 
this example central differencing is used, and Hessian updating is allowed, nonce that tn 

Figure A2.7(d) the Hessian approximation has converged. 

Figure A2.7(e) shows the values of the parameter sensitivity derivatives that were 

calculated by RQSEN. The gradient of the objective function df /dp was calculated by 

1RQCRE ot RQSEN will accept either "YES”, "YE", " 7". "yes , ye.or y for a yes responce 
'’MO" ”N" "no" "n” , null for a no responce. . . . ■ a*** 

2if US CT responded yes then the final aSysis performed using the 

3 Thf* Data file can now be used as an input to RQSEN tor a sensiuvuy dim., 
gradients and Hessian approximation that were found w hen solving e pro 
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— app— , an va.es are » 

pro^U.s^r^pen— 

*** A2.7(e&f))- m «, of - san^p.e oacpa. .nvoWes 

•enni n ati^.eprograi.^ the optimal ^ 0 | S ^^ C pre^tedinFigure 

made. An example of the sensitivity for problem . whe pi 

A2.8. 
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q input name of the file to hold the sensitivity data 

r ts02 

FILEDEF 9 DISK RTS02 DATA Al ( LRECL 800 

INPUT NAME OF THE PROGRAM TO BE RU1 BY RQSEN 

LOAD^RTS02 SRI 4 RQOPT OPTQP RQSEN PARi EN2 PARVAR WRTPT ( CLEAR START 

Execution begins... ********* 

*********************************************** 

* entering rqopt/rqsen preprocessor 
************************************* ********** 


INPUT NAME OF THE FILE WITH INPUT DATA (XXXXXXX) 

OPENING FILE, RTS 02 ON UNIT= 3 FOR INPUT OF PROBLEM DATA 
ARE FUNCT, H, U, V, ACTPT, DELFHG IN FILE RTS02 (YES/NO) 

RQOPT VERSION 2.12 - EXPERMENTAL 2/9/38 


STARTING INFORMATION 


INPUT PARAMETERS : 

NUMBER OF VARAIBL SS = 3 

OPTIMIZATION TO BE PERFORMED ON A NCN-LINEAR OBJECTIVE FUNCTION 
NUMBER OF EQUALITY CONSTRAINTS = 0 

NUMBER OF LINEAR EQ CONSTRAINTS = 0 

NUMBER OF INEQUALITY CONSTRAINTS = 2 

NUMBER OF LINEAR INEQ CONSTRAINTS = 2 

MAXIMUM NUMBER OF ITERATIONS = 100 
LINES PER P/GE = 50 

] TR = 1 

NUMBER OF OUTPUT FILE = 6 


INDEX 


DELTA = 0 . 500000E+00 
R = 0 . 100000E+02 
GAMMA = 0.100000E-04 
EPSILON FOR THE QP = 0.100000E-09 
EPSILON FOR THE GRAD = 0.100000E-03 
DIFFERNCING TYPE = 2 

INITIAL VAULE OF CAPF = O.OOOOOOE+OO 
MINIMUM VALUE OF CAPF = 0.100000E-06 
MINIMUM NORM OF P VECTOR = 0.100000E-03 
SCALING PARAMETER NSC.ALE = 0 

NUMBER OF PARAMETERS NP ARM = 2 

PARAMETER VALUE 

1 1.00000000000000000 

2 3.00000000000000000 


I 


1 

2 

3 


XMIN 

0.000000E+00 

0.000000E+00 

O.OOOOOOE+OO 


XMAX 

0 . 100000E+0 3 
0 . 100000E+C3 
0.100000E+C3 


SCALE 

0 . 100000E+01 
0 . 100000E+01 
0 . 100000E+01 


Figure A2.7(a) A Sample of RQSEN for Test Problem 2. 
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THE HESSIAN APPROXIMATION IN LDL (T) FORM 
ROW 1 0.100000E+01 

ROW 2 0 . OOOOOOE+OO 0.100000E+01 

ROW 3 0. OOOOOOE+OO 0. 000000E+00 0.100000E+01 


ENTERING DRIVER ROUTINE 

ITERATION 

OBJECTIVE FUNCTION = 0.331600C 


PAGE 1 
INDEX 


H (I) = ?0 


V(') 


X(I) 

1 0 . 1100000E+01 

2 0 . 1200000E+01 

3 0 . 1300000E+01 

CHECKING CONVERGENCE THE NORM OF P- 1 
ITERATION 

PAGE 1 OBJECTIVE FUNCTION = 0.255353 

INDEX X(I) = 

1 0. 9059701E+00 

2 0 . 1000000E+01 

3 0 . 1094030E+01 

CHECKING CONVERGENCE THE NORM OF P 0. 
ITERATION 

PAGE 1 OBJECTIVE FUNCTION = 0.255033 

INDEX X(I) H(I) = ?0 v 

1 0 . 1027985E+01 

2 0 . 1000000E+01 

3 0.9720149E+00 

CHECKING CONVERGENCE THE NORM OF P- 0 
ITERATION 

PAGE 1 OBJECTIVE FUNCTION = 0.25500 

INDEX X(I) H (I) = ?0 v 

1 0 . 1000000E+01 

2 0 . 1000000E+01 

3 0 . 1000000E+01 

CHECKING CONVERGENCE THE NORM OF P- 0 
CONVERGENCE ACHIEVED 

FINAL STATISTICS 
CONVERGENCE ACHEIVED 


0 : 0 

OE+02 FUNCTION EVALUATIONS- 1 
CONSTRAINT EVALUATIONS= 14 

g ( i) >= ?o ua) 

0 . 600000E+00 0. 000000E+00 

0 . 170000E+01 0. 000000E+00 


96446237582566496 

1 : 0 

56E+02 FUNCTION EVALUATIONS= 8 
CONSTRAINT EVALUATIONS= 19 
D G (I) >= ?0 U(D 

A-.288658E-14 0. 000000E+00 

A0 . 488060E+00 0. 000000E+00 

345110459418981358 

2 ; 

33E+02 FUNCTION EVALUATIONS= 16 
CONSTRAINT EVALUATIONS= 23 

X) G (I) >= ?0 U(D 

A0. 57731 6E-14 0.586085E+01 

A0 . 244030E+00 0.659360E-01 

395769395140752175E-01 

3 : 2 

J00E+02 FUNCTION EVALUATIONS= 23 
CONSTRAINT EVALUATIONS= 25 
(I) G (I ) >= ?0 U(D 

AO . 288658E-14 0.120000E+02 

A0 . 300000E+00 - . 288658E-14 

. 395476598531791093E-08 


IN 


WITH 


29 FUNCTION EVALUATION 
4 FUNCTION GRADIENTS 
25 CONSTRAINT EVALUATIONS 
2 CONSTRAINT GRADIENT'S 
o TTTERATIONS 

0 . 00000000E+00 BEING THE MAXIMUM CONSTRAINT VIOLATION 

Figure A2.7(b) A sample of RQSEN f :>r test problem 2. 
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THE HESSIAN APPROXIMATION IN LDL(T) 1'ORM 
ROW 1 0.450000E+01 

ROW 2 0 . 444444E+00 0.211111E+01 

ROW 3 0 . 111111E+00 0 . 842105E+00 0. C94737E+01 

THE HESSIAN APPROXIMATION UPPER TRIANGLE 
ROW 1 0 . 450000E+01 0.200000E+01 0 . 50 OOOOE+OO 

ROW 2 0 . 300000E+01 0.200000E+01 

ROW 3 0 . 450000E+01 

PAGE 1 OBJECTIVE FUNCTION = 0 . 2550 J000E+02 FUNCTION EVALUATIONS= 29 

CONSTRAINT EVALUATIONS= 25 
INDEX X(I) H (I) = ?0 V(I) G (I) >= ?0 U(I) 

1 0 . 1000000E+01 A0.288658E-14 0.120000E+02 

2 0 . 1000000E+01 AO . 300000E+00 0.000000E+00 

3 0 . 1000000E+01 

DO YOU WISH TO SAVE THE FINAL DATA (YCS/NO) ? 

y 

DO YOU WISH TO USE THE SAME FILE RTS )2 (YES/NO)? 

n 

INPUT NAME OF THE FILE FOR STORAGE OF DATA (XXXXXXX) 

rts02s 

OPENING FILE, RTS02S ON UNIT= 4 TO STORE PROBLEM DATA 

DO YOU WANT FUNCT, H, U, V, DELFHG, ACTPT WRITTEN TO RTS02S (YES/NO)? 

y 

DO YOU WANT TO PERFORM A SENSITIVITY ANALYSIS (YES/NO) ? 

y 

WELCOME TO RQSEN1 . 0 

A SENSITIVITY ANALYSIS PROGRAM FOR RQOPT 
LAST MODIFIED APRIL 28 1988 

DO YOU WANT TO WRITE THE SOLUTION POINTS TO FILE= 9 (YES/NO)? 

n 


THE DERIVATIVE OF THE OBJECTIVE FUNCTION WITH 
RESPECT TO ALL PARAMETERS WILL BE CALCULATED 

INPUT EPSP FOR THE CALCULATION OF DF / DP 
? 

.0001 


THE DERIVATIVE OF OBJECTIVE FUNCTION W.R.T. P 
-0 . 6999999989E+01 0 . 1200000000E+02 

DO YOU WANT TO STUDY FINITE PERTURBATIONS 

ENTER PARAMETER NUMBER OR (-1 OR CTRL Z) TO CALCULATE GRADS? 
? 

-1 

DO YOU WISH TO FIND PARTIALS OF THE CESIGN VARIABLES AND 
LAGRANGE MULTIPLIERS W.R.T. PARAMETER (NUMBER OR -1 TO END)? 

•? 

1 


Figure A2.7(c) A Sample of RQSEfv for Test Problem 2. 



ENTER EPSP FOR THE GRADIENT CALCULAT] ON? 

? 

.0001 

PERFORMING A SENSITIVITY ANALYSIS FOR PARM ( 1) 

ASSUMING BASE POINT IS STABLE 

BASE POINT VALUE PARM ( 1 ) =0 . 100001 iOOOOOOOOOOE+OI 

PERTURBED VALUE OF PARM ( 1) =0 . 1000: 00000000000E+01 

ENTER THE NUMBER OF ITERATIONS FOR RQOPT? 

? 

2 

************************** entering kqopt ********************** 


ENTERING DRIVER ROUTINE 

ITERATION 0: 0 PARAMETER ( 1) = 0 . 10001000E+01 

PAGE 1 OBJECTIVE FUNCTION = 0 . 255 00500E+02 FUNCTION EVALUATIONS= 1 

CONSTRAINT EVALUATIONS= 26 
INDEX X (I ) H (I) = ?0 V(I) G (I ) >= ?0 U(I) 

1 0 . 1000000E+01 A0.100000E-03 0.120000E+02 

2 0 . 1000000E+01 A0 . 299900E+00 O.OOOOOOE+OO 

3 0 . 1000000E+01 

CHECKING CONVERGENCE THE NORM OF P= 0 . 924648797898831832E-03 

ITERATION 1: 0 PARAMETER ( 1) = 0 . 10001000E+01 

PAGE 1 OBJECTIVE FUNCTION = 0 . 25499301E+02 FUNCTION EVALUATIONS= 14 

CONSTRAINT EVALUATIONS= 28 
INDEX X (I ) H (I) = ?0 V(I) G (I) >= ?0 U(I) 

1 0 . 1000483E+01 A-.133227E-14 0.119995E+02 

2 0 . 9992333E+00 A0 . 299400E+00 O.OOOOOOE+OO 

3 0 . 1000183E+01 

CHECKING CONVERGENCE THE NORM OF P= 0. 673599522571675425E-03 

ITERATION 2: 0 PARAMETER ( 1)= 0 . 10001000E+01 

PAGE 1 OBJECTIVE FUNCTION = 0 . 254 99300E+02 FUNCTION EVALUATIONS= 21 

CONSTRAINT EVALUATIONS= 30 
INDEX X(I) H (I ) = ?0 V(I) G (I) >= ?0 U(I) 

1 0 . 1000208E+01 A-.222045E-14 0.119995E+02 

2 0 . 9997833E+00 A0 . 299400E+00 0.000000E+00 

3 0 . 9999083E+00 

*************************** leaving rqopt ********************** 

THE HESSIAN APPROXIMATION UPPER TRIANGLE 
ROW 1 0 . 500000E+01 0.100000E+01 O.ICOOOOE+Ol 
ROW 2 0 . 500000E+01 0.100000E+01 

ROW 3 0 . 500000E+01 

PERTURBED VALUE OF PARM ( 1 )= 0.999899999999999997 

ENTER THE NUMBER OF ITERATIONS FOR RQOPT? 

1 

Figure A2.7(d) A Sample of RQSEN for Test Problem 2. 
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★ ★★ 


*★★★★★**** 




*** ENTERING RQOPT 


*********** 


*★★★★★***** 


"LSESST'S o "^0 <eIo2 22 

PAGE 1 OBJECTIVE FUNCTION - 0.255 CONSTRAINT EVALUATIONS 

FrtijE, X > G(I) >= ?0 U(l) 

H<1) ' ? ° V<1, S:SSg£S 0 o'.00 0 0 0 0 0 0olt°00 

2 0.1000217E+01 

3 0 . 1000092E+01 0 0 DOOOOOOOOOOOOOOOOE+OO 

CHECKING CONVERGENCE THE NORM OF P ******** 

*“^**=-*-** LEAVING RQ ° Pr ************** 

CENTRAL DIFFERENCE A pp ^ XIAM ^ T ^ N | 9 999 9 c 8 889772823 
DF/DP (PART F + U PART G) = _ 7 ; oo24 99f 3049339675 

- - 1 •<>° 2499 '■ 305129313,, 

D ^ D 2 ^"orT21667Sl , -0.9 1 66 ) 155B + 00 

DV/DP (FINITE DIFFERNECE) 

DU ( 1)/DP( 1) CCE ^ nTFF^ = o’ 00000C 00E+00 
DU ( 2)/DP( 1) (CEOT DIFF) °^0° = . 0 _ 60002999E+01 
ACTIVE CONSTRAINT DG 2 /DP J = . , 0 . 60002999E+01 

PGPP+PXPP*POPP = DG ( 2) /DP ( CHANGE FOR XNCREASE P 

LINEAR ESTIMATE OF WHEN = 0.49998E-01 

G ( 2) ENTERS THE ACTIVE^SET^FOR^D ^ * 0 . l0 499975E+01 

___ TrrrT change for decreased p 

LINEAR ESTIMATE OFWHENACTIVE SET p = _ 0 . 46153E+00 

:i - — - - 

^PERFORMING A PARAME ITER STUDY P™« 

ASSUMING BASE POINT IS STABLE 

BASE POINT VALUE PARM( l,-.l000000j00000000»01 

DF /DP - -0 . 70000000E+01 _ 0 . 216( 7 085E + 01 -0. 91669155EF00 

DX/DP ( 1)= 0-208 x^m vnn THE PARAMITER OR 

S ™ ofSTlO EXIT FORM THIS SUBROUTINES 


.1 


Figure A2 


7(e) A Sample of RQSEN for Tes. Problem 2. 
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THE NEW VALUE OF P ( 1)= 0 . 110000CK E+01 

************************** entering rqopt ********************** 


ENTERING DRIVER ROUTINE 

ITERATION 0: 0 PARAME' 

PAGE 1 OBJECTIVE FUNCTION = 0.248' 

INDEX X (I) H (I) = ?0 

1 0 . 1208319E+01 

2 0 . 7833292E+00 

3 0 . 9083308E+00 

CHECKING CONVERGENCE THE NORM OF P= 
DERATION 


'ER( 1)= 0 . 11000000E+01 

•3909E+02 FUNCTION EVALUATIONS= 30 
CONSTRAINT EVALUATIONS= 60 
V(I) G (I) >= ?0 U(I) 

A0.208111E-01 0 . 115333E+02 

A-.300030E+00 0.000000E+00 


PAGE 1 


0.219176271755135058 

1: 0 PARAMETER ( 1) = 0 . 11000000E+01 

OBJECTIVE FUNCTION = 0 . 247 42404E+02 FUNCTION EVALUATIONS= 37 

CONSTRAINT EVALUATIONS= 62 
V(I) G (I) >= ?0 U(I) 

A0.415223E-13 0.104876E+02 

A0 . 235367E-13 0.542741E+00 


H (I ) = ?0 


INDEX X(I) 

1 0 . 1045910E+01 

2 0 . 7944067E+00 

3 0 . 1055092E+01 

CHECKING CONVERGENCE THE NORM OF P : 
CONVERGENCE ACHIEVED 
*★***★*******■★***★★★*★**★★* LEAVING 
ENTER THE PERTURBATION FOR THE P ARAMS 
ENTER 0.0 OR NULL TO EXIT FORM THIS 
? 

0.0 

DO YOU WANT TO STUDY FINITE PERTURBA' 
ENTER PARAMETER NUMBER OR (-1 OR CTR 

? 

-1 

DO YOU WISH TO FIND PARTIALS OF THE 
LAGRANGE MULTIPLIERS W.R.T. PARAMETE 

•? 

-1 

Ready; T=2.81/4.12 11:23:12 $2.46 


O.OOOOOOOOOOOOOOOOOOE+OO 

:ter or 

SUBROUTINE? 


:ions 

j Z) TO CALCULATE GRADS? 


)ESIGN VARIABLES AND 
* (NUMBER OR -1 TO END) ? 


Figure A2.7(f) A Sample of RQSEN for Test Problem 2. 
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